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\Ve  review  and  unify  those  applications  and  techniques  associated  with  Feynman’s  theory  of  path 
integrals  which  have  been  found  relevant  for  sound  propagation  in  the  ocean.  After  giving  an  intro- 
ductory discussion  of  functional  integrals  in  general  and  Feynman  path  integrals  in  particular,  we 
derive  several  path  integral  representations  for  the  solutions  to  the  two-  and  three-dimensional  parabolic 
equations.  The  analogies  which  exist  be'ween  sound  propagation,  the  nonrelativistic  quantum  mechan- 
ics of  a point  particle,  and  Brownian  motion  are  con^ered.  Next  we  use  the  path  integral  to  derive  , , 

several  methods  of  approximation  including  perturbation  theory,  the  Rytov  approximation,  ray  acoustics,  ^ ^ 
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and  straight-line  geometric  optics.  The  formalism  is  then  applied  to  the  problem  of  developing  algo- 
rithms for  numerically  solving  the  parabolic  equation.  After  developing  path  integral  representations 
for  the  solution  to  the  Helmholtz  equation,  we  give  an  extensive  discussion  of  the  parabolic  approx- 
imation. 
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AN  INTRODUCTION  TO  THE  APPLICATION 
OF  FEYNMAN  PATH  INTEGRALS  TO  SOUND  PROPAGATION 

IN  THE  OCEAN 

1.  INTRODUCTION 

In  an  increasing  number  of  situations,  one  is  interested  in  determining  the  features  of 
acoustic  propagation  in  a general,  range-dependent,  ocean  environment.  By  way  of  illustration, 
during  a recent  meeting  of  the  Acoustical  Society  of  America  HI  55%  of  the  talks  devoted  to 
underwater  propagation  dealt  with  phenomena  which  cannot  be  adequately  modeled  in  terms 
of  an  ocean  medium  possessing  range-independent  characteristics. 

Of  the  techniques  which  are  available  for  modeling  propagation  in  a general  environment, 
the  parabolic-equation  technique  of  Tappert  and  Hardin  [2]  is  perhaps  the  most  promising.  It 
has  a wider  region  of  validity  than  ray-tracing  techniques,  it  takes  into  account  mode  coupling 
to  ail  orders,  and  it  can  be  implemented  numerically  using  relatively  simple  algorithms. 

The  technique  assumes  that  the  solution  to  the  Helmholtz  equation  for  the  pressure  field 
is  not  so  much  different  from  the  solution  to  a simpler  equation  called  the  parabolic  equation. 

Because  the  Helmholtz  equation  is  an  elliptic  differential  equation,  one  must  simultaneously 
solve  for  the  field  at  all  points  in  order  to  obtain  the  field  at  any  particular  point.  This  requires 
considerable  computation  for  a range-dependent  medium.  On  the  other  hand,  the  structure  of 
the  parabolic  equation  allows  one  to  obtain  its  solution  at  a particular  range  knowing  only  its 
solution  at  shorter  ranges.  The  equation  can  therefore  be  solved  by  marching  in  the  sense  that 
the  field  at  range  r can  be  obtained  by  propagating  the  field  at  a somewhat  shorter  range 
r — Ar  according  to  a version  of  Huygens’  principle.  The  computational  realization  of  this  con- 
cept leads  to  the  above-mentioned  algorithms,  the  best  known  of  which  is  the  split-step 
Fourier  algorithm  |2,3|. 

Most  of  the  past  interest  in  the  parabolic  equation  has  centered  about  this  computational  ' 

technique.  The  equation  was  almost  always  discussed  in  the  context  of  some  scheme  for  ob- 
taining a numerical  solution.  It  was  realized  only  recently  that  the  equation  is  also  a useful 
starting  point  for  analytic  studies  of  sound  propagation. 

Of  the  analytic  tools  available  for  study  of  the  parabolic  equation,  Feynman's  theory  of 
path  integrals  |4|  seems  particularly  powerful.  It  is  very  intuitive.  One  expresses  the  full  wave 
solution  to  the  parabolic  equation  in  terms  of  the  quantities  associated  with  ray  acoustics.  The 
split-step  Fourier  algorithm  follows  directly  in  a few  lines  from  the  discrete  version  of  the  path 
integral.  In  fact  the  most  efficient  way  of  solving  the  parabolic  equation  is  by  evaluating  the 
path  integral.  The  split-step  algorithm  is  just  one  of  several  techniques  for  doing  this.  The  pri- 
mary defect  of  Tappert  and  Hardin's  technique  is  that  a realistic  ocean  bottom  cannot  be 
modeled.  This  is  a fault  of  the  algorithm  rather  than  the  basic  parabolic  approximation. 

Feynman’s  theory  points  the  way  toward  the  development  of  an  algorithm  for  computing  the 
pressure  field  above  and  in  a physical  ocean  bottom.  The  parabolic  approximation  itself  is  most 
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naturally  discussed  using  the  path-integral  lormalisni  /inally,  Dashen  (5|  has  recently 
developed  a unified  theory  of  sound  propagation  through  a random  ocean  by  using  path  in- 
tegrals. 

Ihc  purpose  of  this  report  is  to  dispel  some  of  the  mystery  which  surrounds  path  in- 
tegrals and  their  application  to  sound  propagation.  We  have  attempted  to  present  the  introduc- 
tory material  which  anyone  interested  in  using  path  integrals  as  a research  tool  would  need  or 
want  to  know.  In  essence  the  report  consists  of  a collection  of  examples  and  applications. 
I here  does  not  seem  to  be  any  way  of  becoming  comfortable  with  path  integrals  without  work- 
ing through  many  examples.  Although  the  report  is  rather  mathematicai,  containing  over  500 
equations,  the  mathematics  which  is  used  is  pedestrian.  The  whole  subject  is  treated  at  an  ele- 
mentary level,  l-'or  those  interested  in  a more  sophisticated  discussion  of  Feynman’s  theory, 
several  excellent  review  articles  are  available.  Among  these  we  have  found  the  works  by  Mon- 
troll |6|,  Ciel’fand  and  Yaglom  [71,  Brush  18|,  Tarski  |9|,  Kravtsov  |10|,  Brittin  and  Chappell 
111),  I'radkin  (12),  Berry  and  Mount  (13),  Klyatskin  [14],  Keller  and  McLaughlin  (151,  and 
Koeling  and  Malfliet  (I6|  particularly  useful.  The  standard  reference  is  the  monograph  by 
Feynman  and  Hibbs  |17].  In  addition  to  a general  exposition  of  the  subject  this  work  contains 
a host  of  useful  ct^mputational  techniques. 

This  report  is  organized  into  nine  sections  and  two  appendixes.  Section  2 contains  some 
comments  concerning  functional  integrals  in  general  and  Feynman  path  integrals  in  particular 
After  f'eynman  path  integral  representations  arc  derived  for  'he  two-  and  three-dimensional 
parabolic  equations  in  Section  3,  several  alternative  representations  are  developed  In  Section  4. 
In  Section  5 the  analogies  are  pointed  out  which  exist  between  sound  propagation,  quantum 
mechanics,  and  Brownian  motion.  In  Section  6 path  integrals  are  used  to  derive  several 
methods  of  approximatif)n.  Br)undary  conditions  are  discussed  in  Section  7,  and  the  utility  of 
path  integrals  for  developing  algorithms  is  emphasized.  The  solution  to  the  Helmholtz  equation 
IS  written  as  a path  integral  in  Section  8,  and  several  applications  are  considered.  Section  9 is  a 
discussion  of  the  parabolic  approximation.  This  report  does  not  specifically  consider  the  applica- 
tion of  Feynman’s  theory  to  the  problem  of  sound  propagation  through  a random  medium,  be- 
cause the  subject  has  been  treated  with  thoroughness  by  Dashen  (5).  However,  so  that  some 
previous  work  can  be  discussed,  acoustic  fluctuations  are  briefly  considered  in  Section  9.  Ap- 
pendix A contains  the  derivation  of  an  intergral  representation  for  the  ii-dimensional  Green’s 
function  in  an  infinite,  homogeneous  medium.  This  representation  is  used  throught  the  report. 
In  Appendix  B the  stationary-phase  approximation  of  an  integral  is  recalled. 

The  various  sections  of  the  report  are  not  intended  to  be  read  separately.  Each  section 
uses  previously  derived  results. 


2.  FUNCTIONAL  INTEGRALS 
2.1  General  Comments 

A functional  integral  is  a generalization  of  an  ordinary  A-dimensional  integral.  To  obtain 
the  value  of  the  A/-dimensional  integral 

“ i"  *^^1  dZf^J{2\ 2f^),  (2.1) 
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one  evaluates  the  integrand  J over  the  range  of  the  N variables  Z| , , , These  values  of  J 

are  then  added  according  to  the  rules  of  calculus  to  determine  1/^.  A functional  integral 
corresponds  to  the  continuum  limit  of  If^.  The  index  i of  z,  becomes  a continuous  variable,  say 
r. 

(Z|,Z2 Zyy)  — |z(/-|  ),  z(r2  h ....  z(/-^))  — zf/-).  /-continuous.  (2,2a) 

The  integral  over  the  /V-tuple  lz| . Z;^.}  becomes  an  integral  over  the  function  zfr): 

/ t/z,  • • • t/z^  - / Z)[z(r)l,  (2.2b) 

And  the  integrand  7 becomes  a functional,  that  is  a function  of  a function: 

y(z| . ...,  z^)  — ' y(z  (/■) ).  (2.2c) 

One  writes  for  the  functional  integral 

/=/ Z)[z(/-)1  y(z  (/■)).  (2.3) 

Operationally  it  means  7 is  evaluated  for  all  permitted  functio  (r)  and  then  these  contribu- 
tions are  added  according  to  the  rules  of  functional  calculus  to  gi  ■ of  the  integral. 

There  are  three  important  points  concerning  functional  integi,.  irst  point  is  that 

J(z(r))  is  not  simply  a function  of  the  parameter  r,  it  cannot  be  determ  5y  picking  a value 
of  r,  evaluating  z at  that  point,  and  then  determining  7 for  that  value  ol  z.  The  functional  7 can 
be  determined  only  if  the  complete  function  z(r)  is  specified,  not  its  value  at  any  one  point. 
An  example  of  a functional  which  illustrates  this  point  is  the  area  A(z(r))  under  a curve  z(r) 
as  shown  in  Fig.  1.  Obviously  A has  a value  which  cannot  be  determined  unless  the  complete 
function  z is  given.  Knowing  z at  some  particular  value  of  r is  not  enough.  It  should  be  clear 
from  this  example  that  J(z(r))  is  the  same  functional  as  7(z  (s) ) and 

/ D(z(r)l  7(z(r))  =/  D[z  (.s)  1 7(z  (v) ). 

Our  notation,  though  traditional,  is  somewhat  misleading  in  this  regard.  In  the  mathematics 
literature  one  often  finds  the  preferable  notation  7(z(  • ) ) and  J D|z(  • ) 1 7(z  ( • ) ). 
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Fig.  I — The  area  under  a curve  as  an  example  of  a functional; 
A(z(r))  = f dr  z(r). 
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The  second  point  is  that  one  must  know  the  class  of  functions  involved  to  determine  the 
value  of  a functional  integral.  This  corresponds  to  the  fact  that  the  /V-dimensional  integral  can- 
not be  evaluated  unless  the  limits  of  integration  are  given.  The  particular  specification  of  the 
class  of  functions  depends  on  the  problem,  but  some  specification  must  always  be  given  As  an 
example,  one  might  specify  that  the  integration  is  over  all  continuous  and,  of  course,  single- 
valued functions  z(r)  defined  on  the  interval  n K r ^ (i  such  that  ;(«)  =C|,  z(fi}  = ^2’ 
and  A ^ z{r) 

The  third  point  concerns  how  the  various  contributions  are  added  to  form  the  integral. 
Integral  calculus  is  just  the  study  of  the  rules  by  which  this  addition  is  accomplished.  Different 
sets  of  rules  correspond  to  different  types  of  integration,  e.g.,  Riemann-Stieltjes,  Lebesque,  etc. 
These  rules  determine  how  the  contributions  to  the  integral  are  to  be  weighted;  they  determine 
the  measure.  The  measure  must  always  be  specified.  An  expression  such  as  hq.  (2.3)  is  mean- 
ingless unless  the  measure  is  explicitly  or  implicity  indicated.  Without  it  one  would  have  no 
idea  how  to  go  about  determining  the  integral's  value. 


2.2  Feynman  Path  Integrals 

In  1948  Feyn.man  I-*'  sublished  a formulation  of  nonrelativistic  quantum  mechanics  con- 
ceptually distinct  from  the  )rmulations  developed  during  the  mid-1920’s  by  Schrodinger  and 
Heisenberg.  In  Feynman’s  formulation  the  wave  equation  and  operator  calculus  are  replaced 
by  a type  of  functional  integration  called  path  integration. 


Specifically,  consider  the  one  dimensional  motion  of  a particle  of  mass  m moving  under 
the  influence  of  a position-  and  t'  re-dependent  force  FU,  z).  If  this  force  is  derivable  from  a 
potential. 


FU.z)  “ -d,y(t.z). 


(2.4) 


then  the  Lagrangian  for  the  system  is 

L m 1-^1  - F(f,  z) 

2 \di 


kinetic  energy  — potential  energy. 


(2.5) 


According  to  Hamilton’s  Principle  [18]  the  particle  will  move  from  z^  at  time  to  z^  at  time  tf, 
in  such  a way  that  the  action 


A = J di  Liz,  i) 


(2.6) 


is  an  extremum.  This  requirement  leads  to  the  differential  equation  for  the  path  followed  by 
the  particle  iNewton’s  law): 


m 


dr  2 


-f),y(i.z)  = FU.z). 


(2.7) 


^herezlf,,)  -z„andz(r*) 


If  now  we  consider  the  particle  to  be  described  by  quantum-mechanical  rather  than  clas- 
sical dynamics,  its  motion  is  characterized  by  a probability  amplitude  'i’  il^,  z^]//,,  Zf,).  The 
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modulus  of  'I'  gives  (he  probability  (hat  the  particle  will  travel  from  at  time  to  Z/,  at  time 
Feynman  showed  that 

'F  = X (2.8) 

paths 


where  A is  given  by  Eq.  (2.6),  H is  Planck's  constant  divided  by  In,  and  the  sum  is  over  all  the 
possible  paths,  not  just  the  classical  one(s),  connecting  (/„,  z^)  to  (//,,  Z/,).  Equation  (2.8)  may 
be  rewritten  as 


:a\'h.  exp 


exp 

1 

— ni 

a 

2 

(tl 

y(i,z(i)) 


(2.9) 


The  complete  quantum-mechanical  behavior  of  the  particle  is  summarized  by  Eq.  (2.9). 
Although  no  mention  need  be  made  of  the  dil.erential  equation  satisifed  by  'I'  or  of  the  com- 
mutation relations  satisfied  by  the  observables,  hcynman  demonstrated  that  his  formulation  is 
mathematically  equivalent  ’o  the  traditional  formulations. 


Because  of  its  close  connection  to  classical  physics.  Feynman’s  theory  is  perhaps  the  most 
intuitive  approach  to  quantum  mechanics.  It  has  become  increasingly  important  in  theoretical 
physics.  For  example,  since  1965  an  average  of  at  least  one  paper  every  nine  days  has  been 
published  in  which  some  aspect  of  Feynman’s  theory  is  considered.  Applications  have  been 
found  in  quantum  field  theory  17,9,12,19-22],  statistical  mechanics  [6-8,11,12,17,23-27],  the 
development  of  asymptotic  expansions  for  solutions  of  differential  equations  1 1 3,1 5,28-29],  the 
development  of  variational  techniques  ]l7,30j,  the  formulation  of  semiclassical  approximations 
to  scattering  processes  [16,31-33],  the  determination  of  bound  states  [34],  the  quantization  of 
field  equations  [35-39],  the  study  of  disordered  systems  [40-47],  and  the  somewhat  related 
problem  of  wave  propagation  through  a random  medium  [5,10,14,48,491  f he  theory  has  gen- 
erated interest  among  mathematicians  in  the  problem  of  rigorously  defining  path  integrals 
]5()-55]  and  in  the  relationship  between  functional  integrals  and  differential  equations  [7,56-59]. 
Moreover  work  has  been  devoted  to  approximation  schemes  and  the  numerical  evalution  of 
path  integrals  [8,60-69].  This  incomplete  list  gives  some  idea  of  the  impact  the  theory  has  had 
on  physics  and  mathematics. 


2.3  The  Question  of  Rigor 

Feynman  path  integrals  have  never  been  given  a rigorous  mathematical  definition  The 
basic  problem  is  that  the  functional  integral  of  Eq.  (2.9)  does  not  have  a countably  additive 
measure  [50].  This  statement  simply  means  that,  with  the  usual  methods  of  measure  theory, 
path  integrals  are  senseless.  Although  there  exist  more-or-less  satisfactory  definitions  of  path 
integrals  [50-55],  they  apply  only  to  certain  forms  of  the  potential  Ff/,  z). 

How  should  an  acoustician  view  this  lack  of  mathematical  rigor?  He  should  probably  ig- 
nore it  and  adopt  (he  point  of  view  that  (he  path  integral  "if  manipulated  in  a purely  formal 
style  without  any  regard  for  rigorous  justification,  gives  all  the  right  answers"  [70].  After  all,  if 
one  dw'.v  find  an  inconsistency  which  is  solely  the  result  of  the  use  of  path  integrals,  he  has 
made  a discovery  of  profound  importance. 
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3.  THK  KKYNMAN  PATIMNTWiRAL  SOIAJTION  TO  THK  FARABOIAC  KQLATION 
3.1.  Some  Preliminaries 


I hroughoul  ilrs  ..urvey  we  assume  one  is  interested  in  determining  the  pressure  field 
[>(x,  I)  at  some  point  x = ix,  y,  z)  in  the  ocean  and  at  some  time  i due  to  a point  source  locat- 
ed at  x^  = z,).  One  traditionally  assumes  obeys  the  wave  equation 


- 6 (x  - X,)  f (i). 


(3.1) 


where  ~ The  function  /(i)  defines  the  spectral  characteristics  of  the 

source,  (f  or  a monochromatic  source,  ,/(r)  =acos((u„/  + </().)  The  temporal  variability  of 
the  sound  speed  clx,  t)  is  usually  gradual  enough  so  that  Kq.  (3.1)  can  be  replaced  by  the 
llelmholt/  equation 


2 

-E 

OJ 

c(x,  /) 

- fi  '■'’'  (\  - X,  ), 


(3.2) 


where 


pix,  1)  = 

j d<„  Ki<»)p„,ix,  t)  e 

(3.3a) 

fil) 

= f tJtu  f;(o))c 

(3.3b) 

In  addition  to  satisfying  the  Helmholtz  equation,  (x,  t)  will  satisfy  the  appropriate  boundary 
conditions  at  the  lop  and  bottom  of  the  ocean  and  will  represent  outgoing  radiation.  With  re- 
gard to  the  solution  of  Eq.  (3.2),  both  w and  I enter  as  free  parameters,  i.e.,  the  equation  is 
solved  as  if  they  were  constant.  We  shall  therefore  suppress  the  dependence  of /)  on  these  vari- 
ables and  write  pix)  for  the  complex  pressure  i). 


fn  the  parabolic  approximation  one  assumes 

p(x)  = const  — ^ z),  (3.4) 

« 1/2 

where  R - (x  - -E  (y  - y^)^  and  -=  vj/c„,  with  c„  being  some  reference  sound 
speed.  The  reduced  field  i//  satisfies  the  parabolic  equation 


-lik„b,*\i(r,z)  = lf)j  ■{-  kl  \n^  (r,  z)  - 1 11  i/<  (r,  z), 
where  «(r,  z)  -=  t„/(  (r,  z),  with 

r , , . / 


cir.  z)  = c 


X,  + — (x  -x,).y,  ' 

A A 


(3.5) 

(3.6) 


The  variable  r is  the  horizontal  path  length  in  the  direction  of  propagation  (It  is  commonly 
said  that  the  parabolic  approximation  is  related  .somehow  to  cylindrical  symmetry.  Strictly 
speaking  this  is  not  correct.  What  the  parabolic  approximation  implies  is  that  the  acoustic  en- 
ergy is  restricted  to  the  vertical  slice  of  the  ocean  containing  the  source  and  receiver.  Since 
only  this  slice  is  important,  the  three-dimensional  problem  reduces  to  a problem  involving  only 
two  coordinate  variables:  rand  z.) 
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The  solution  to  Eq.  (3.5)  satisfies  the  same  boundary  conditions  at  the  top  and  bottom  of 
the  ocean  satisfied  by  p(x)  and  an  initial  condition.  This  initial  condition  is  usually  specified 
by  requiring  that  i/i  equal  some  depth-dependent  function  at  a small  horizontal  distance  « from 
the  source.  That  is, 

ill  (f,  z)  = h{z).  (3.7) 

The  precise  form  of  h and  the  value  of  f depend  on  the  constant  coefficient  in  Eq.  (3.4),  the 
acoustic  properties  of  the  ocean  in  the  region  of  the  source,  and  one’s  inclination.  For  any 
choice  of  initial  field,  however,  iji  can  be  written  in  terms  of  a Green’s  function  4': 

ill  (r.  z)  = J dzQ  (r,  z\€,  Zq)  hiz^),  r ^ €,  (3.8) 

where 

( +2/)tydf  + dj  +k^[n^(r,z)  — 1 ])  'F  (r,  z]  r',  z')  =0  (3.9a) 

and,  for  any  value  of  r', 

'F(r',  z|r',  z')  = S (z  - z').  (3.9b) 

This  Green’s  function  propagates  the  acoustic  energy  from  a source  at  r',  z'  to  the  point  r,  z. 
We  always  take  r > r'in  Eq.  (3.9a). 


If  the  effects  of  the  variation  of  the  sound  speed  in  cross  range  are  considered  important, 
one  could  use  the  three-dimensional  version  of  the  parabolic  equation.  Taking  the  direction  of 
propagation  to  be  along  the  positive  x axis,  we  have 


p(x)  = const  V ' ' ^ix  — X,,  f)}. 

(3.10) 

where  p 

is  a two-dimensional  vector,  p - (y.  z). 

ihlr.f'i}  - dp„  'P  (r,  pI*.  p„))i(p„), 

(3.11) 

1 + 2/*,,  + k !.i’(r,  p)  - 111  'Kr,  p|r.  p)  = 0, 

(3.12a) 

and 

'P  (r'.  p|  r',  p ')  — fi*^’(p  — p ), 

(3.12b) 

with  = by  + 9?  and  n(r.  (i)  - c„/(  (x,  -t-  r,  y,  z,  r). 

3.2  The  Huygens-Fresnel  Principle 


Of  the  many  methods  available  for  developing  the  Feynman  path-integral  solution  to  the 
parabolic  equation,  perhaps  the  easiest  to  understand  is  the  one  based  on  the  Huygens-Fresnel 
Principle  (4).  We  will  illustrate  this  method  by  considering  the  solution  to  Eq.  (3.5)  in  an  un- 
bounded medium. 

The  Huygens-Fresnel  principle  gives  a technique  for  constructing  the  field  at  r in  terms 
of  the  field  at  a shorter  range  r',  where  r — r'  = Aris  taken  to  be  small.  Every  element  of  the 
field  at  r' is  considered  the  source  of  secondary  wavelets  which  coherently  sum  to  produce  the 
field  at  r.  These  wavelets  are  solutions  to  the  governing  wave  equation  with  constant  index  of 
refrzction.  That  is,  the  rays  associated  with  the  wavelets  follow  straight  paths.  Let  z,  be  a 
depth  variable  at  range  rand  let  z,  be  a depth  variable  at  r' as  illustrated  in  Fig.  2.  (It  might 
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Fig.  2 - The  Huygens-Fresnel  Principle 


be  more  customary  to  simply  use  z and  z'  for  z^  and  z,  respectively.  We  chose  this  notation 
because  it  will  suggest  a functional  integration.)  According  to  the  Huygens-Fresnel  principle 
(71)  the  contribution  to  i/<  (r,  z^)  from  an  element  at  r',  z,.  is 

rfi//(r,z,)  = T)  (r',  z^)  wCr,  z^  ; r',  z,)i/i  (/•',  z,.  )rfz,’,  (3.13) 

where  w(r,z^  -,  r'.z^  ) is  the  field  associated  with  the  wavelet.  The  presence  of  the 
phenomenological  factor  tj  is  characteristic  of  the  Huygens-Fresnel  theory;  it  is  usually  defined 
within  the  context  of  the  Kirchhoff  formulation.  It  will  be  determined  below  by  a consistency 
condition. 


Now  w(r,  z, ; /■',  z,  ) satisfies  Eqs.  (3.9)  with  n^(r,  z,)  taken  to  be  constant.  Writing 

i(r—r')k„  , ) 

w - exp  (rt^(r,z,) -1  ij  K,  (3.14a) 

we  have 

(+2/k<,d,  + )w  -0  (3.14b) 

and 


- 6(Zf  - Zr')- 

Eqs.  (3.14)  are  easily  solved  by  Fourier  transforming  with  respect  to  z^.  We  obtain 


(3.14c) 


w{r,  z,  ; r',  z,  ) 


27r/(r  — /•') 


iko  ir  - r')  z,  - z. 


+ fl2(r,z,)  - 
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Integrating  Eq.  (3.13)  over  all  the  elements  gives 

L '(2  I 

r(  I 

<l>(r.z,)  ” - ; — J dz, 

2ir;(r  — r ) •' -oo  ' I 


ikjr-r')  7,  -7,. 


+ 7,)  - 1 


X T)  (r',  7,.)i/«  (r’,  z^'). 


The  factor  y)  is  determined  by  the  condition  that  Eq.  (3.16)  reduce  to  an  identity  as  r—r'.  In 
this  limit 

i/<  (r,  7,)  — J"  rf7,.  T)  (/■,  7^.)  8 (7^  - 7,  ) 0 (r.  7,.) 

— 7)  (r,  7,)  0 (/■,  7,), 

implying 

T,  - 1.  (3.17) 

(Equation  (3.17)  is  also  obtained  from  the  parabolic-equation  analog  of  the  Kirchhoff  theory.) 

It  is  convenient  to  define 

U{r,z)  = -~[n'^{r.z)  -1).  (3.18) 

\f  nir,  z)  =1  — K (r,  7),  where  x is  a small  perturbation,  then  U{.r,  7)  = k (r,  z)  . 

In  summary,  the  field  at  r is  given  in  terms  of  the  field  at  r'  by 


C 

= -T — I dz,&X\> 


/7c„Ar  — 


2 Ar 


— (/(r.  7^)  I 0 (r',  7,  ),  (3.19) 


where  \r  r — r'. 


By  carrying  out  the  indicated  differentiations,  one  can  show 
{-2ik„df  - - *0  lr»^(r.  7,)  -lll0(r,  7,) 

— - — —likgb^n^ir.Zf)  —b^n^ir.z^) 

— — ^ — [dj  n^(r,  7,)  1 ^ 0(r,  7,)  - Ib^  n^  (r,  z^)b^  <i)  (r,  z,) 

If  Ar  is  chosen  small  enough  so  that  the  gradients  of  n^(r,  7^)  can  be  ignored,  the  right-hand 
side  of  Eq.  (3.20)  is  insignificant  and  Eq.  (3.19)  is  a good  approximation  to  the  solution  of  Eq. 
(3.S).  In  other  words,  if  Ar  is  picked  small  enough  so  that  the  rays  associated  with  each 
wavelet  do  not  get  a chance  to  significantly  bend  as  a result  of  the  variation  of  the  index  of  re- 
fraction, Eq.  (3.19)  will  give  a valid  approximation  to  the  field  at  rand  7,. 

By  an  iterative  procedure  we  will  now  obtain  the  field  at  range  r •=  R and  depth  7 due  to 
a source  at  r — 0 with  the  spatial  distribution  at  r — * given  by  Eq.  (3.7).  We  partition  the 
plane  of  propagation  into  N strips  of  width  Ar  by  lines  at  r " rj , r2,  ....  ryy_|  as  indicated  in 
Fig.  3.  Just  as  before,  a depth  coordinate  z,  is  associated  with  the  boundary  r ~ r,.  The  nota- 
tion is  simplified  if  we  set 


(3.20) 
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Fig.  3 — A single  Feynman  path  z(r).  The  summation  over  all  the  paths 
amounts  to  integrating  over  Zf^,  Zf2,  ■■■,  2ryv_i 

The  Huygens-Fresnel  principle  will  now  be  applied  to  each  strip.  At  r — /q  we  have, 
from  Eq.  (3.7), 

('■q,  (3.22a) 

At  the  end  of  the  first  strip 

\ k I 

|2 


\ 2ir/Ar 
X expl/ACgAr 


z..  - z,. 


Ar 


(3.22b) 


At  the  end  of  the  second  strip 


k 


X exp\ikg^r 
I 2/2 


1 

2,  “ 1 

2| 

Ar  1 

27r/Ar| 


X exp 


/ dz^  h{z^  ) 
I l^r, 


ikgl^r  £ 
/-I 


Ar 


~ ^r,  ) 


(3.22c) 
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At  the  end  of  N strips 


k 

(R.  z)  - / dz,^  dz,^  • • • h(z,^  ) 


^\\ 


X exp  ik„£kr  T — 

" ,h12  Ar 


- f/(r,,  2,  ) . (3-22d) 


Because  of  the  finite  size  of  Ar,  Eq.  (3.22d)  is  only  an  approximation.  The  exact  solution 
is  obtained  as  Ar  — * 0.  Since  AAr  — /?  — t is  fixed,  as  Ar  — » 0 — » <»  and  the  AZ-dimensional 

integral  in  Eq.  (3.22d)  becomes  infinite  dimensional.  This  infinite-dimensional  integral  is  inter- 
preted as  a functional  or  path  integral.  The  points  (/-q,  z,^),  (r^,  z,^)  ....  (r^,  z,J  are  con- 
sidered to  lie  on  a continuous  curve  z(r)  such  that  2,  - z{r,),  giving  z(e)  - zCrg)  - zg  and 
z(r/v)  - z(R)  - 2 (Fig.  3).  The  A - 1 intermediate  integrations  over  2, , ...,  2,  become  an 
integration  over  z(r): 


fd^r,  dzr,.,  - S D[z(r)]. 


This  equation  defines  the  measure  associated  with  the  path  integral.  The  expression 
(z,^  ~ z,  )/A^r  — (2,  — 2,  ^ )/(r,  — _l  ) is  the  discrete  approximation  to  the  derivative  of 

2(r),  and  Ar]^  (...)  is  the  discrete  approximation  to  an  integral  over  r.  Therefore,  as 


Ar  — * 0, 


N , 2,-2, 

. ^ 1 '•  ''-I 

i T ^ 

2 Ar 


R r 2 

- i/(r,  z,  ) - r rfr  - Uir.zir))  . (3.24) 

' ' J % 1 dr 


Collecting  expressions  gives 


f f 2 

4>(R,z)  ^ J dzQ  hizg)  f Dlz(r)]exp  '^0  J* -U{r,z(r))  , (3.25) 

where  the  integration  is  over  all  continuous  functions  2(r)  such  that  2(e)  - Zq  and 
z{R)  — 2.  Comparing  Eq.  (3.25)  with  Eq.  (3.8)  gives 

\ 1 1 

'F  (/?,  2|«,  Zq  ) - J*  Z)[z(r)  1 exp  '*0  J"  -i/(r,  z(r))  (3.26) 

for  the  Green’s  function.  For  an  arbitrary  coordinate  dependence 

'F  (r.  2I  r',  2')  - J*  Z)l2(s) ) exp  J* j"^  —U(s,z(s))  , (3.27) 

where  z(s)  is  a continuous  curve  defined  on  the  interval  r' < 5 < rand  satisfying  the  end- 
point conditions  2 (r')  — z'andz(r)  — z. 


II 
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Equations  (3.25),  (3.26),  and  (3.27)  are  exact;  we  have  "only"  assumed  the  path  integral 
exists. 

3.3  The  Composition  Law 

In  this  subsection  we  will  outline  a second  technique  for  constructing  the  path  integral. 
We  will  work  with  the  Green’s  function  defined  by  Eqs.  (3.9).  It  satisfies  the  composition  law : 

(r,  z\r',  z‘)  - J*  t/Zj  9' (r,  zj  s,  Zj ) 't' (s,  Zjl  r',  z'),  (3.28) 

where 

s^  r‘.  (3.29) 

In  the  context  of  the  study  of  Markov  processes,  Eq.  (3.28)  is  sometimes  referred  to  as  the 
Smoluchowski-Kolmogorov  relation. 

To  show  that  the  composition  law  is  valid,  we  must  show  that  the  right-hand  side  of  Eq. 
(3.28)  satisfies  Eqs.  (3.9).  Equation  (3.9a)  follows  easily  enough  by  operating  on  both  sides  of 
Eq.  (3.28)  with  2/)c„d,  -I-  b}  + [n^ir,  z)  — 1 ) and  noting  that  this  operator  can  be  moved 
inside  the  integral  over  z^.  To  show  Eq.  (3.9b)  holds,  one  must  take  the  limit  r—  /-'without 
violating  Eq.  (3.29).  This  can  be  done  in  two  ways:  r— • s followed  by  s — r'or  s — r' followed 
by  r — s.  In  either  case  Eq.  (3.9b)  will  be  satisfied.  This  two-step  limiting  process  is  artificial 
and  can  be  avoided  by  doing  the  integral  over  Zj  for  r - sand  s - r' small  but  finite  and  then 
taking  the  limit  r—  r'.  Below  we  will  obtain  an  expression  for  for  small  horizontal  separa- 
tions which  permits  this  procedure  to  be  carried  out.  The  result  is  again  Eq.  (3.9b). 

It  is  rather  interesting  that  the  right-hand  side  of  the  composition  law  is  independent  of  s. 
This  characteristic  is  a result  of  the  falloff  of  'I'  at  infinity.  To  see  this,  we  differentiate  the 
right-hand  side  of  Eq.  (3.28)  with  respect  to  sand  use  the  equation  reciprocal  to  Eq.  (3.9a), 

\+lik^b,.  -bl  -klW(r'.z')  - 111 '!'(/■,  zlr'.z ')  -0,  (3.30) 

to  obtain 

dj  J*  dz^  S'  (r,  z|  s,  Zj ) 'I'  (s,  z,]  r',  z') 

“ liT  S Id/,  A s.  ) '»'  (i.  2, 1 r‘,  z') 

— Mr  (r,  z|  s,  Zj  )d/  Mr  (s,  z,|  /■',  z')  ] 

- Mr  (r.  z\  s,  Zj ) Mr  (s,  zJ  r'.  z') 

- Mr  (/-,  z|  s,  z, ) bj^  Mr  (s.  Zj  I /■',  z ') )/'  (3.3 1 ) 

-0, 

provided  the  Green’s  function  falls  off  at  infinity  (or  provided  the  Wronskian-like  quantity  in 
the  last  pair  of  brackets  is  independent  of  z^). 

We  now  iterate  the  composition  law 

Mr  (r,  z|r',  z')  - J dz^  . . .dz^_^  0 Z/|S( -I • z, -i ).  (3.32) 

(-1 
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where 

Sy  = r,  Zy  = z,  Sq  = r',  Zq  = z'.  (3.33) 

We  pick  s,  (/  = 1 , N - 1 ) so  that  s - s,  _)  = As,  giving  N^s  = r - r'.  Consider  one  of 
the  factors  in  Eq.  (3.32).  Since  ^ (s,  z,  js,  , z,  _, ) = (s,  z,|  s,  - As,  z,  _, ) is  a 8-function 

in  z,  — z,  for  As  = 0,  it  will  be  "almost  like  a 8-function"  for  As  small.  This  means  that  only 
a small  region  around  z,  = z,  _j  contributes  to  the  integral  in  Eq.  (3.32).  Therefore,  in  calcu- 
lating 'I' (s;  z,ls;_l,  z,_i ),  we  need  only  consider  the  variation  of  the  index  of  refraction 
within  a small  region  about  the  point  s,,  z,.  In  particular,  by  appropriate  choice  of  As,  this  re- 
gion can  be  made  so  small  that  n(r,  z)  will  be  essentially  constant  within  it  and 

( k.  \ik.  (fz.  -z,_.  P 

'•'(s,,  z,ts^  -1-Z/-1  ) 


ko 

1/2 

-^As 

2,  -2/-1 

27r/As 

exp 

2 

A5 

-1-  /I^(s,,  Z,  ) - 1 


,(3.34) 


which  solves  Eq.  (3.9)  for  a constant  index  of  refraction.  Substituting  Eq.  (3.34)  into  Eq.  (3.32) 
and  using  Eq.  (3.18),  we  obtain 

I W2 

/ </z,  rfeyv-i 


'I'  (r,  z|  r',  z') 


X exp 


N 

iko  As  £ 
/-I 


-2,-1 


As 


- (/(s,.  z,  ) 


(3.35) 


for  As  small.  Just  as  in  the  treatment  following  Eq.  (3.22d),  Eq.  (3.35)  reduces  to  a functional 
integral  as  As  — • 0: 

(zq,  Z) zjyl  — z(s),  (3.36a) 

where,  from  Eq.  (3.33),  r'  < s < r,  z(r')  = z',  and  z(r)  = z, 

N/2 

j dz^  dZy.i  - / f)[z(s)l 

|2 


Ini^s 


As^ 
, -1 


1 

2,  2/ 

2 

As 

- t/(s,,  z,) 


^ \ ds 

1 

dz 

j 

r 

2 

ds 

giving  Eq.  (3.27)  again: 

r 

'Ir  (r,  z\  r'.  z')  - J*  D[z  (s)  1 exp  ikg  f ds 

r‘ 

3.4  The  Three-Dimensional  Parabolic  Equation 


- i/(s,  z(s)) 


— f/(s,  z(s) ) 


(3.36b) 

(3.36c) 


(3.27) 


The  derivation  of  the  path-integral  solution  to  Eqs.  (3.12)  follows  step  by  step  the  deriva- 
tion for  the  two-dimensional  problem.  One  obtains 


'V  {r,  p\r',  p')  - J"  Dip  (s)l  exp 


r 

iko  J 


- 


2 ds 


U(s,  p(s)) 


(3.37) 
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The  integration  is  over  a two-dimensional  vector  function  p(4)  such  that  p (r  ) 
P (r)  = p.  The  measure  is  defined  by  the  rule 


27r;As 


N 

f 


IrriLs 

/ Z)[p(s)l. 


N 

S ^>'1  dz\ 


i>'  and 


(3,38) 


4.  ALTERNATIVE  REPRESENTATIONS 


In  the  previous  section  we  derived  a path-integral  solution  of  the  form 


(4.1) 


where  /4  is  a functional  of  the  path  z(.s)  and  its  derivative.  This  representation  is  usually  called 
the  coordinate  representation.  For  many  problems  it  is  not  particularly  convenient.  In  this  sec- 
tion we  develop  two  alternative  representations. 

4.1.  The  Wavenumber  Representation 

In  the  wavenumber  or  momentum  representation  the  number  of  degrees  of  freedom  is 
doubled  by  introducing  an  additional  integration  over  a function  k(s).  Physically  k(s) 
represents  the  vertical  wavenumber  associated  with  the  path  z(s). 

We  first  rewrite  Eq.  (3.34)  in  the  form 


/till, 


'l'(s,.Z,|s,_,.Z,_,  ) - 

- y z,)  I 

(We  have  defined  A,  to  be  dimensionless ) Substituting  Eq.  (4.2)  into  Eq.  (3.32)  yields 
'P(r,  z|r',  z ) ^dz^  ..  dz^_^^dk^  dfc^^cxpj 


(4.2) 


N 

X Z 

I -I 


-Z/-I 


As 


(*„As 


In  the  continuum  limit 


As  ^ 


Izq-  ^/vl  “* 

1*1  , ...,  */vl  -•  *(s), 


Zi  - 2 


/ -I 


As 
- S ds 


- T ■2/  ^ 


k{s)  ~ --^*2(s)  -l/(s,  2(s)) 
ds  2 


(4.3) 

(4.4a) 

(4.4b) 

{4.4c) 
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giving 


I* 

^ S <ih- S ■ ■ dk^-^  j D[z(s)]Dik(s)] 

'V  (r,  z\  r',  z')  - J"  £)lz(s)]Dl/c(s)l  exp|//c„ 

X 'jds\k(s)-^  - t/(s,z(s))|j. 


The  path  k(s)  is  unconstrained,  and,  as  before,  z(r')  — z'and  z(r)  — z.  From  a comparison 
of  Eqs.  (4.4d)  and  (3.36c)  we  notice  a different  measure  is  being  used  for  the  functional  in- 
tegral in  the  wavenumber  representation. 

In  three  dimensions  the  wavenumber  representation  is 


S' (r,  p|r'.  p ) - J Dlp(s)lDU(s)l  exp  /*„ 


X j <ys|if(s)  ■ - U(.s,  p(s))j|, 


where 


2N 

J rfp,  ..  dkf^-^  J Dlfi(s)]Dlk(s)]. 


4.2  The  Velocity  Representation 

Many  times  it  is  desirable  to  transform  a path  integral  by  employing  a change  of  variable; 
z(s)  — * z'(s)  “■  some  function  of  z(s). 

We  shall  introduce  this  technique  by  transforming  Eq.  (3.27)  into  a representation  in  which 
the  integration  is  over  the  "velocity" 

, . dz 
v(s)  “ — . 
ds 

We  first  observe  that  Eq.  (3.35)  can  be  written  in  the  form 

'i' (r,  z|  r',  z')  - / </z,  ...  dz^6(z  - Z;v)  exp  |/)c„As 


We  now  introduce  a change  of  variable 


z,  -Z(-i 


, / - 1 


D.R  PALMER 


or,  since  Zq  “ z', 

Z]  “=  z'  + AiV] , 

Zj  “2'  + AsCv)  + V2 ), 


Zjy  “ z'  + As(v|  + V2  + ...  + P/v). 

The  Jacobian  of  the  transformation  is  simply 

</z,  ...  r/z/y  ”=  (^s)’^dv^  ...  dvj^. 
Substituting  this  into  Eq.  (4.8)  gives 


(4.10) 

(4.11) 


4'  Ir.  z|  r’.  z ) 


\ Itti 
N I 


m 


s 


dvi 


dvf^  6 z - z'  - As  ^ vJ  exp  {//c^^As 


/-I 


^ 2 
— y*  — 

2 ' 


f/  s,,  z'  + A.V  ^ Vj 


/-I 


(4.12) 


In  the  continuum  limit 

r 

'I' (r,  z|  r',  z')  “ J £)lv(s)l8  z - z'  - 5 ds  v(s) 


exp  liky 


X j* 


y v2(i)  - U 


S,  z'  + J*  c/.v'  v(4') 


(4.13) 


Unlike  the  previous  examples,  the  integration  is  over  a function  v(s)  which  is  not  constrained 
by  end-point  conditions.  The  measure  is  determined  by  the  correspondence 

W2 

J"  dv,  ..  dv/^—  ^ Dlv(s)l.  (4.14) 


2ni 


For  the  three-dimensional  problem 

4/ (r,  pjr',  p')  = J D(v(4)1  8 - p'  - J rfs  v(s) j exp 


-y  v^(s)  - U 


s.  p ' -F 


where 


xfds 

r' 

[■^fl  J ‘^i’l  J j)lv(s)l. 


,1 

J*  ds'  Us') 

r' 


ik„ 


(4.15) 


(4.16) 


The  velocity  representation  occurs  frequently  in  the  Russian  literature. 

4.3.  Normalization  Conventions 


We  have  defined  the  measure  by  indicating  that  the  functional  integral  is  the  limit  of  a 
multidimensional  integral  as  the  dimensionality  approaches  infinity  and  then  specifying  the 
normalization  of  this  multidimensional  integral.  That  is,  Eq.  (3.36b)  indicates  how  Eq.  (3.27)  is 
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to  be  evaluated.  This  procedure  is  usually  not  followed  in  research  articles.  Instead  the  normali- 
zation is  implicity  specified  by  stating,  for  example,  that  the  integral  in  Ecj.  (3.27)  is  defined 
such  that 


J Dlz(s)  1 exp 


2 

2 ^ 

^ t 

ds 

1/2 

ikjr  - r ) ^ 

2 

2Tr/(r  — r') 

exp 

2 r - r 

(4.17) 


We  want  to  show  in  this  subsection  that  this  procedure  is  entirely  equivalent  to  the  one  we 
have  been  using.  First,  Eq.  (4.17)  is  certainly  valid  if  we  assume  the  conventions  given  by  Eq. 
(3.36).  By  comparing  Eqs.  (4.17)  and  (3.27)  it  is  clear  the  left-hand  side  of  Eq.  (4.17)  is  the 
solution  to  Eqs.  (3.9)  with  n — 1.  This  solution  can  be  trivially  obtained  directly  from  Eqs. 
(3.9)  and  is  equal  to  the  right-hand  side  of  Eq.  (4.17). 

We  now  want  to  show  that  Eq.  (4.17)  does  indeed  give  the  normalization  indicated  in  Eq. 
(3.36b)  by  explicity  evaluating  the  path  integral  in  Eq.  (4.17)  The  left-hand  side  of  Eq.  (4.17) 
is  first  written  as 


\ . oo 

where,  as  before,  zq  ~ z',  Zyy 


ik„  ^ 

J 

rfz/v_|  exp 

z,  and  A.$  = (r  — r')/N.  We  want  to  show 
I m 


AiN)  ' 

Following  the  steps  which  led  to  Eq.  (4.12), 


k„N 


2ni(r  ~ r') 


(4.18) 


(4.19) 


ik„ 

N 

J 

rfz;y_,  exp 

0 

2As 

L - z,  _, ) ^ 

/ -1 

\ 

//c„ 

s 

J 

. . t/v/v  ^ 

z - z'  - As  ^ V, 

/ “I 

exp 

0 

2As 

I 

- ’ n / exp 


2 As 


2 _ 


i&,skv, 


- 

ko 

dk 

(2v 

— exp 

ik(z  ■ 

-z')  - 

A'As  .2 

\ 2ni(r 

--■>1 

m 

1 1 

1/2 

exp 

ik„ 

27ri  (r 

- r') 

2(r  - r'] 

z) 


In  writing  Eq.  (4.20)  we  have  used  the  important  formula 


1/2 


,0^/An 


(4.20) 


(4  21) 


Comparison  of  Eqs.  (4.20),  (4.18),  and  (4.17)  gives  (4.19),  the  desired  result. 
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We  have  followed  a standard  scheme  m evaluating  the  path  integral;  a change  of  variable 
is  introduced  to  diagonali/e  the  exponent,  and  the  integrations  are  then  carried  out  using  the 
formulas  for  Gaussian  integrals. 


5.  SOME  ANALOGIES 

5.1  Quantum  Mechanics 

The  quantum-mechanical  behavior  of  a particle  of  mass  m moving  along  the  z axis  in  a 
potential  yi.t,z)  is  determined  by  a wave  function  if  (/.z)  which  satisfies  Schrodinger’s  equation 

-2i  ~ d,>l>  U.z)  = d}  y(i.  z)  0.  z).  (5.1) 

where  il  is  Planck’s  constant  divided  by  In.  If  the  wave  function  is  known  at  some  initial  time 
Iq,  it  will  be  given  at  a later  time  by 

<li(i.z)  = J dzQ  K(i,  z\iq,  Zq)  \l)  (Iq,  Zq),  (5.2) 

where  the  Green’s  function  kernel  (Feynman  propagator)  K obeys  the  Schrodinger’s  equation 
and  the  initial  condition 

AC  (Iq.  z| /q,  Zq  ) =8(z  — Zq).  (5.3) 

If  one  compares  Eqs.  (5.1),  (5.2),  and  (5.3)  with  Eqs.  (3.5),  (3.8),  and  (3.9),  one  sees  there  exists 
an  analogy  between  the  nonrelativistic  quantum-mechanical  dynamics  of  a point  particle  and 
acoustic  propagation  in  the  ocean.  The  correspondence  can  be  made  precise  by  introducing 
some  reference  speed  c„  which  may  be  though  of  as  the  average  speed  of  the  particle.  Its 
value  is  unimportant,  since  it  will  cancel  in  any  quantum-mechanical  expression.  It  is  intro- 
duced simply  to  make  the  dimensions  work  out.  In  Table  1 we  have  constructed  a dictionary 
for  passing  from  quantum  mechanics  to  underwater  sound  propagation.  The  meaning  of 
several  items  in  the  table  will  become  clear  when  we  consider  the  ray-acoustics  approximation 
to  the  parabolic  equation. 

It  may  be  a mistake  to  attach  great  importance  to  this  analogy.  One  cannot  solve  the 
problems  of  sound  propagation  in  the  ocean  by  simply  borrowing  quantum-mechanical  results. 
Much  work  has  been  devoted  to  quantum  mechanics,  however,  and  techniques  and  specific  cal- 
culations which  are  now  being  discussed  by  acousticians  first  appeared  in  a quantum- 
mechanical  context  and  have  been  available  for  years. 

5.2  Brownian  Motion 

If  a particle  constrained  to  move  along  the  z axis  undergoes  Brownian  motion,  the  posi- 
tion of  the  particle  as  a function  of  time  is  described  by  a probability  density  >l>(i,z)  which 
satisfies  the  diffusion  equation 

a,  •/«  <5.4) 


The  diffusion  coefficient  22  expresses  the  physical  characteristics  of  the  particle  and  medium  via 
Einstein’s  relation 
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Table  1:  A Quantum  Mechanics- Acoustic  Propagation  Dictionary 


DR  PALMBR 


t i 

Flz(s)]  = exp  “ “7/  "X  £>^  U(s)]. 


The  probability  that  the  particle  followed  name  path  is  equal  to  the  integral  of  Eq.  (5.7)  over  all 
paths  and  must  be  given  by  Eq.  (5.5): 

+x.- 

/ dz'  J Plzis)]  = Pit.  z)A2,  (5.8) 

z paths 

where  the  sum  is  over  all  paths  such  that  z(0)  — z^  and  z(i)  = z'.  Equation  (5.8)  implies 


L f 

dz 

2 

1 

1/2 

(/  -/q) 

^ 

2 

-s-° 

r-i 

ds 

2-n Dh  — /q  ) 

6 Xp 

2D 

' -'0 

where  the  integration  is  over  all  paths  z (s)  such  that  z(Iq)  = Zq  and  z (r)  = z. 

Consider  now  the  calculation  of  the  average  value  of  a functional  E(z(s))  of  the  path 
followed  by  the  particle.  This  average  is  given  by  the  functional  integral 

r D^lz(s)]exp  S F(z(s)) 

2D  ds 


< f(z(v))  > 


/D,.(z(x))exp 


1*1 


If  we  take 


P(z(s))  “ exp r t/.s  f'fs,  z(s) ), 


we  have 


J_ 

dz 

D 

ds 

< F(z(s) ) > - iV  r /)  [z(s)  1 exp  — ?:  f ds  — ~ - K(.v.  z(s) ) , (5.12) 

D 2 ds 


where  N\sa.  normalization  constant 

/V-'  - //)Jz(.s)lexp  • 


The  functional  integrals  Eqs.  (5.9),  (5.10),  and  (5.12)  are  called  Wiener  integrals  (721.  All 
the  statistical  aspects  of  Brownian  motion  can  be  expressed  in  terms  of  these  integrals. 


From  the  form  of  Eq.  (5.12)  we  see  the  analogy  between  Wiener  integrals  and  Feynman 
path  integrals.  A Feynman  integral  can  be  (formally)  converted  into  a Wiener  integral  by  ana- 
lytic continuation  in  the  mass,  in  the  time,  or  in  ft.  Moreover,  since  they  can  be  rigorously 
defined,  Wiener  integrals  have  often  been  the  starting  point  of  studies  attempting  to  rigorously 
define  Feynman  path  integrals. 
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The  brevity  of  the  preceding  discussion  really  does  not  do  justice  to  the  application  of 
Wiener  integrals  to  the  theory  of  Brownian  motion.  References  7,  8,  and  17  contain  richer, 
more  satisfying  discussions. 


6.  METHODS  OF  APPROXIMATION 


In  virtually  every  realistic  situation  the  parabolic  equation  cannot  be  solved  analytically 
and  one  must,  therefore,  resort  to  some  method  of  approximation.  The  methods  which  are 
used  in  acoustics  generally  fall  within  two  categories:  asymptotic  methods,  such  as  ray  acoustics, 
and  perturbation  methods,  such  as  the  Born  approximation  or  the  Rytov  approximation.  In  this 
section  we  will  use  path  integrals  to  develop  two  methods  belonging  to  each  category. 


6.1  Straight-Line  Geometric  Optics 


If  acoustic  energy  is  propagated  from  (0,  Zj)  to  (R,  z),  Eq.  (3.27)  becomes 


'P  (R,  z|0,  z^)  - J"  Z)(z(r)  1 exp 


ik. 


A 


i- 

2 dr 


- U(r.z(r)) 


where  z(0)  -z^andzC/?)  — z.  The  integral  is  normalized  so  that 
J*  D\z(r)  \ exp 


fdr 

dz 

2 

*0 

1/2 

ikgR 

z -V 

2 

2 V' 

dr 

2mR 

exp 

2 

R 

(6.1a) 


(6.1b) 


In  straight-line  geometric  optics  the  assumption  is  made  that  a ray  picture  is  valid  and  refrac- 
tion is  unimportant.  In  such  a situation  the  acoustic  phase  depends  only  on  the  ray  path  which 
follows  a straight  line  from  source  to  receiver.  The  equation  for  this  path  is 


Vo.  ('■>  “V  + (z  - ZjJ-  ^6.2) 

This  suggests  introducing  a change  of  path  variable 

z(r)  —Zgg(r)  + z'(r)  (6.3) 

and  integrate  over  z'(r)  rather  than  z(r).  Since  z(0)  = z^„  (0)  = Zj  and  z(R) 
z^g  (R)  ” z,  the  function  z' satisfies  the  endpoint  conditions 

z'(0)  -z'(R)  -0.  (6.4) 


Moreover,  the  Jacobian  of  the  transformation  is  unity.  To  see  this,  we  return  to  the  discrete 
version  of  the  path  integral.  Equation  (6.3)  becomes 

-z.,  (z  -zj  +z;.  (6.5) 

where  / “ 0 N.  Therefore 

J r/z,  rfz^_,  - J r/z,  <6.6) 

since  Eq.  (6.5)  corresponds  to  simply  adding  constants  to  the  variables  z,.  Equations  (6.1)  be- 
come 
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M' (/?,  z|0.  ) -=  exp 

- u 


,k„R 

2 

Ik  (5. 

— f dr 

1 1 

dz' 

2 

R 

1 u \z  \»/j  exp 

2 V' 

[2 

dr 

r,  Zj  + (z  — z,)  + z'(r) 

A 


(6.7a) 


and 


J Dlz'(r)  1 exp 

r dr 

dz' 

2 

I 

z'{R)  -z'(0) 

dr 

ImRl 

exp  1 2 

R 

2niR 


1/2 


(6.7b) 


So  far  we  have  made  no  approximations.  We  now  assume  that  in  the  path  integration  of  Eq. 
(6.7a)  the  dependence  of  f/on  z'can  be  ignored.  That  part  of  the  integrand  depending  on  U 
can  then  be  moved  outside  the  integral  sign,  and  the  integral  over  z'  is  just  the  normalization 
condition  Eq.  (6.7b).  We  obtain 


'K/?,  zjO,  z,)  ='«/„(/?,  z|0,  Zj)  exp 


K 

- { dr  U r.  z,  + -^  (z  - z, ) 


(6.8) 


where  is  the  solution  to  the  parabolic  equation  with  n(r,z)  — 1 (i.e.  with  U - 0): 


z|0,  z,) 


J^o_ 

1/2 

ik„R 

2 

z -z* 

liriR 

exp 

2 

R 

(6.9) 


The  straight-line  geometric-optics  approximation  is  important  for  the  study  of  electromag- 
netic propagation  through  the  atmosphere  but  it  is  of  little  value  for  underwater  sound  propa- 
gation (at  moderate  to  long  ranges),  because  the  sound  channel  always  introduces  significant 
refraction.  However  this  example  does  illustrate  an  important  methodological  point.  Many 
times  it  is  not  necessary  to  actually  evaluate  a path  integral  in  order  to  obtain  a desired  result. 
In  this  example  the  normalization  condition  determined  the  only  path  integral  needed. 


6.2  Ray  Acoustics 

The  path  integral  can  be  viewed  as  a coherent  sum  of  elementary  ray  solutions 
exp  \ik„A(z(r)\  corresponding  to  various  paths  z(r).  These  ray  paths  do  not  satisfy  any  type 
of  eikonal  or  ray  equation  however.  They  are  completely  arbitrary  with  the  exception  they 
must  start  at  the  source  and  end  at  the  receiver.  Consider  what  happens  when  k„  becomes 
large.  Most  of  the  paths  will  contribute  little  to  the  sum,  because  exp  [ik„A  (z(r)l  will  oscillate 
rapidly  and  be  canceled  by  the  equally  rapidly  oscillating  contributions  from  other_ nearby  paths. 
The  only  exceptions  will  be  for  those  paths  z(r)  which,  when  perturbed  (z  — • z -h  z),  do  not 
produce  any  significant  change  in  k„A.  The  difference  between  A(z  (r)  -t-  r(r))  and  A(z  (r)) 
will  never  be  zero  for  all  paths  z , but  it  will  be  sufficiently  small  if  it  vanishes  to  first  order  in  z 
and  if  terms  higher  than  second  order  are  small.  The  situation  is  entirely  analogous  to  that 
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encountered  when  one  evaluates  an  integral  by  the  method  of  stationary  phase.  The  only 
difference  is  that  one  deals  with  an  infinite-dimensional  integral  rather  than  a finite- 
dimensional integral. 


Turning  now  to  specifics,  we  introduce  in  Eqs.  (6.1)  the  change  of  variable 

z(r)  = z*(r)  + z‘(r) 

with 

z'iO)  = z'(R)  =0 

and  expand  the  exponent  through  second  order  in  z': 

(A(.  z|0,  Zj ) = exp  \ik„A{z*(r))]  J*  Dlz'fr)  1 exp 

R 

i 


(6.10a) 

(6.10b) 


R 

y ( dr 


- z'9,t/(r.  z‘)  -b  ~ 
dr  dr  ~ 2 


dz‘ 


dr 


(z')^d}U{r.z') 

2 


where 


first-order  cancellation  requires 


AU-(rn  -t,,,.,-, 

res 


(6.11) 


(6.12) 


or,  integrating  by  parts  and  using  Eq.  (6.10b), 

d^z' 
dr^ 


n 

/ drz'ir) 


+ d.U(r.  z') 


0. 


The  only  way  Eq.  (6.13)  can  be  satisfied  for  all  z'  is  for  z*  to  satisfy  the  equation 

d^z' 


dr^ 


= -Q.Uir.z*). 


(6.13) 


(6.14) 


Equation  (6.14),  together  with  the  endpoint  conditons  z*(0)  = z,  and  z*{R)  = z determine 
the  ray  path  z*.  Equation  (6.14)  corresponds  to  Newton’s  second  law  F = ma.  We  now  rewrite 
Eq.  (6.11): 


where 


'I' (/?,  z|0,  z^ ) = exp  (z*(/-)  1 T, 

T =«  J”  Dlz'(r)  ] exp  | r/r||-^|  — (z')  ^Mir) 


— 1* 

dz' 

2 

dr 

with 


Mir)  = d^Uir.z'ir)). 


(6.15a) 


(6.15b) 


(6.16) 
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Again,  the  path  integral  is  normalized  so  that 


f Dlz  (r)] 


exp 


— ] dr\ 

dz‘ 

2 

ko 

2 •! 

dr 

lit  iR 

1/2 


(6.17) 


Before  proceeding  with  the  evaluation  of  T,  let  us  indicate  how  the  straight-line 
geometric-optics  approximation  is  recovered.  We  take  U(r.  z)  to  be  such  a slowly  varying  func- 
tion of  z that  d;U(.r,  z^  in  Eq.  (6.14)  and  Mir)  in  Eq.  (6.15b)  can  be  ignored.  The  solution 
to  Eq.  (6.14)  is  then  Zy„  ir)  of  Eq.  (6.2)  and  we  immediately  obtain  Eq.  (6.8). 


The  integral  Tcan  be  written  in  the  discrete  form 

T 


where  (dropping  the  prime  on  the  path) 

t m 

f 


lim 

/V  '^op 

ArA  -R 


(6.18) 


Ini^t 


dz 


\-l 


X exp|- 
ik,Ar 


z?M,  + Zt  M^  -I- 


2 ■^^2 '"2  +z^_,A/yv-i 

with  Ml  = Mii\r)  = i/(/Ar,  z*(/A/-)).  Collecting  terms  in  the  exponent  gives 

. 1 

“P  (1^7  I, 

The  matrix  iA^j)  has  components 

A,,,  - 2 - (Ar)2  M,. 

and 

0,  1/  -7|  > 2. 

We  now  introduce  the  subdeterminants 
Do  -1 


(6.19) 


(6.20) 


and 

That  is 


£)*  - det  (A  1),  i.j  - 1 k. 


D, 


(6.21) 


(6.22) 


D,  - 


'N-\ 


^11' 

■^11 

-4,2 

-421 

-422 

" del  A. 
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Qn  ~ 

(6.28) 

i 

* 

SO  that 

K 

27r/0^ 

(6.29) 

i 

\ 

From  Eqs.  (6.21)  and  (6.23)  we  have  the  recursion  relation 

Qk-\-\  ~ Qk  Ok  ~ Ok-i  ^ ^ 

Ar  Ar  * * 

(6.30a) 

for  /c  — 1 N — \ . Moreover 

1 

o 

(6.30b) 

and 

1 

Q\  ~ Oo  j 

Ar 

(6.30c) 

i 

i 

f 
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As  A/-  — * 0, 

IC^o C^/vl  - 0(r).  0 ^ R. 

and  Hqs.  (6.30)  become,  respectively, 

+ M(r)Q(r)  - 0, 

Jr^ 

Q(0)  -0, 

and 


(6.31a) 

(6.31b) 


jJ2 


(Jr 


(0)  - 1. 


(6.31c) 


We  therefore  obtain  the  result 


J li\zir)  1 exp 

\ dr 

dz 

2 

•*2  ( t.\  k4  ( r\ 

1/2 

2 

dr 

— z Kr) M \r) 

liTiQiR) 

(6.32) 


for  z(())  z(R)  —0,  where  Q(r)  satisfies  Eqs.  (6.31).  Equation  (6.32)  is  valid  provided 
Q(R)  ^ 0,  i.e.,  provided  (/4,^)  is  a positive  definite  matrix.  This  example  illustrates  the  close 
connection  between  Feynman  path  integrals  and  boundary-value  problems  of  the  Sturm- 
LiouvilJc  type  J57J. 


We  must  now  solve  Eqs.  (6.31)  for  M(r)  - U(r.  z*(r)).Tht  ray  equation,  Eq.  (6.14), 
may  be  replaced  by  the  equations 


r 

z*(r)  — z,  + J”  ds  pis) 


(6.33a) 


and 


p(r)  “ p(0) 


-|r/.v 


B,U(s.z*(s)). 


where 


p(r)  — 


dz*(r) 


dr 


(6.33b) 

(6.33c) 


Physically  pir)  - tan  dir).  = dir)  where  dir)  is  the  local  angle  the  ray  path  z*ir)  makes 
with  respect  to  the  horizontal.  Consider  the  quantity  dz*ir)ldpiB).  Differentiating  Eq.  (6.33a) 
gives 


dz*ir)  r , dpis) 


dpiO) 


- j ds 


dp  iO ) ' 

But  from  Eq.  (6.33b)  we  have 

-1  - / <^»'6/(/(s-,z*(s'))  4^- 

dpiO)  i ‘ dpiO) 


(6.34a) 


dpiO) 

Substituting  Eq.  (6.34b)  into  (6.34a)  gives  the  integral  equation 

dz*ir)  C . C ^,..5.,/..  dz*is') 


(6.34b) 


dpiQ) 


(6.3W 
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II 

I 


k 


If  we  differentiate  Eq.  (6.34c)  with  respect  to  r,  we  obtain 

Idz'ir) 


and 


dr^ 


dr 

dz*(r) 


- 1 


dp{Q) 


dpm 

b}U(r.  z'(r)) 


- Mir) 


dp  iO) 


dpiO)  ■ 


From  Eqs.  (6.34c)  and  (6.34d)  we  have 


dz*ir) 


and 


dr 


dp  iO) 


dz*(r) 


r-0 


Comparing  Eqs.  (6.34e),  (6.34f),  and  (6.34g)  with  Eqs.  (6.31),  we  see 

dz*(r) 


dp  iO) 
4g) 
Qir) 


-■  1. 


r-0 


dpiO) 

Combining  Eqs.  (6.15),  (6.32),  and  (6.35),  we  finally  find 

i 

'!'(/?,  z|az,) 


InidzIdpiO) 


exp  \ik^Aiz*{r))\, 


(6.34d) 


(6.34e) 


(6.34f) 


(6.34g) 


(6.35) 


(6.36) 


where  p(0)  is  defined  by  Eq.  (6.33c),  A by  Eq.  (6.12),  and  the  ray  path  z*by  Eq.  (6.14). 


If  the  matrix  (A^)  of  Eqs.  (6.20)  and  (6.21)  is  not  positive  definite,  det  A — Qi/f)  — 
dz/dpiO)  —0.  The  endpoint  z*(R)  — z of  the  ray  path  is  insensitive  to  the  initial  grazing  an- 
gle, which  means  the  ray  path  is  degenerate,  and  a caustic  is  present.  This  leads  to  an  infinite 
amplitude,  and  the  ray  acoustics  approximation  breaks  down. 


We  can  go  further  Returning  to  Eqs.  (6.27)  and  (6.28),  we  have 


1 


J_ 

dzidp  (0 ) /v‘-~  I Qj 


“ lim 


Qs- 


lim 


rfz*(r,  )/dpi0) 


Qn 


dz*(r^  )/dpi0) 


dz*ir-)  )/dpi0) 


rfz*(z/v_|  )/dp(0) 


dz*ir^)/dpi0) 


- ( - 1 )'  lim 

N —OB 


where  / is  the  number  of  the  factors 


I 


dz*^r^  )/dpi0) 


dz*ir,_i  )/dpi0) 
dz*iri)/dpiO) 


dz*ir^_^  )/dpi0) 
dz*ir^)ldpiO) 


(6.37) 


(6.38) 


which  are  negative.  If  this  ratio  is  negative,  </z*(/-(_l  )/</p(0)  and  rfz*(r^ )/i/p(0)  are  opposite  in 
sign,  implying  that  somewhere  in  the  range  to  r,,  dz*/dpi0)  —0,  (a  caustic  was  encoun- 
tered). Therefore, as  )V— 


I 

1 

I 
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1 

1/2  . 

— tnl 

1 

dz/dp(0) 

” exp  2 

dz/dp  (0 ) 

1/2 


(6.39) 


where  / is  the  number  of  caustics  encountered  as  one  propagates  from  0 to  R.  Consider  now 
Eq.  (6.12)  for  /4(z*(r)); 

2 

- z,*) I 


A iz*{r) ) — lim  Ar  T 

A 


Ar-*0 


/-I 


Z,  -z,-i 


Remembering  that  zg  “ Zj,  we  have 
a/((z*(r)) 


dz. 


— - lim 


Ar 


Z|  - 2n 


Af 


-/»(0), 


so  that 


1 


_ dp{0)  _ 

dz/Jp{0)  dz 

Combining  Eqs.  (6.36),  (6.39),  and  (6.41)  gives 
'!'(/?,  z|0.  z,) 


d^A  (z*(r)) 
dzSz, 


ko 

d^A  (z*) 

1/2 

2iri 

dzdz, 

exp 

ik„A(z*)  - 


ml 


(6.40) 


(6.41) 


(6.42) 


where  again  exp  (—ivl/2)  represents  the  cumulative  effect  of  phase  jump  of  it/2  on  passing 
through  each  caustic. 

We  have  assumed  until  now  that  a single  solution  z* exists  to  the  ray  equation,  if  several 
distinct  rays  contribute,  'I'  is  given  by  a sum  of  terms  of  the  form  of  (6.42): 


'('(/?,  z|0,z,)  - 1 


2ir/ 


d^A  (z/) 
dzdzj 


1/2 


exp 


ikoA  (Zj  ) ^ 


(6.43) 


(The  integer  Ij  depends  on  the  particular  ray  path.)  For  two  solutions  z,*  and  z’  to  the 
ray  equation  to  be  considered  distinct  rays,  the  phase  difference  kjAiz^  ) — A(zj  )|  must  be 
greater  than  some  fraction  of  a cycle,  say  n/2: 


\A  (z/)  — a {zj)\  > (acoustic  wavelength)/4 . 


(6.44) 


We  have  spent  considerable  space  on  ray  acoustics  because  it  is  probably  the  most  impor- 
tant approximation  and  it  gave  us  the  opportunity  to  evaluate  an  important  path  integral  (Eq. 
(6.32)). 

6.3  Standard  Perturbation  Theory 

We  start  our  development  of  perturbation  methods  by  expanding  that  part  of  the  ex- 
ponential in  Eq.  (6.1a)  which  depends  on  If. 

R 

- ikg  drUir,  z(r) ) 

- Z ^ J 1 (;(r,,z(r, ))  U(r„.z(r„)). 

» -0  "•  i)  0 
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By  introducing  the  Fourier  transform 

' r e''^'U{r.k). 

•'  \27T  ) 

(he  right-hand  side  of  Eq.  (6.45)  becomes 

( - '<(„)"  r /•  dk,  dk„ 


n -0  "•  0 


. C 

(27r)  (277  ) 

R 

X exp  ik^  J*  dr  z(r)Mir), 


(/(r,.k,  ) (/(^.k„) 


where  Mir)  is  a sum  of  8 functions: 

Mir) 


^ k,8  (r  - r, ). 


Substituting  this  expansion  into  Eq.  (6.1a)  gives 

, Z <-'*«)"  r r 

'I'  (/?,  2 0.  r ) “=  T I dr^  dr„  \ — 

' „•%  n\  i'  "Jiln)  (277 ) 

X J" /)|z (r)  1 exp  /k„  J"  </7  — ^ + zir)Mir) 


Uir..k,  ) . . Uir„.  k„) 


iMir)  depends  on  r^,  r„,  k^,  ....  k„.)  Although  we  will  later  evaluate  this  path  integral 

directly,  it  is  convenient  at  this  stage  to  observe  Eq.  (6.49)  has  the  structure 


“ Z -7  jf  ‘^'■i  jf  ' J- 

n-O  "•  0 0 


(6.50a) 


where /is  symmetric  under  permutations  of  the  r/s.  Since  there  are  n!  of  these  permutations, 
Eq.  (6.50a)  may  be  written  as 


R 

'1 

1*4*.-, 

1 

" Z J J -I  J dr,  /ir,.  r„).  (6.50b) 

ft  ■<)  0 ^ 0 

The  result  of  rewriting  the  expression  for  'F  in  the  form  of  Eq.  (6.50b)  is  that  the  variables 
7| , ....  r„  which  appear  in  Mir)  become  ordered. 

0 ^ 7,  ^ 72  < < 7„_,  < 7„  < «.  (6.51) 

With  this  ordering  the  path  integral  in  Eq.  (6.49)  can  be  expressed  in  a particularly  convenient 
form.  Writing  7-4/  iR,  z|0,  z, ) for  this  path  integral,  we  first  note  it  is  the  solution  to  the  para- 
bolic equation  with  Uir.  z)  — — zMir)  and  hence  obeys  the  composition  law.  Suppose  we 
iterate  the  composition  law  and  pick  as  intermediate  range  points  the  variables  r^ 7^; 

T/^fiR.^\0.^,)  - J dZi  dz„T„iR.z\r„,z„) 

^ TMir„.z„\r„.^.z„.\)  T^ir^,  z^\Q.  z,).  (6.52) 

(Remember  the  composition  law  by  itself  is  exact.)  We  now  write  a path-integral  solution  for 
the  factors  in  Eq.  (6.52): 
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'k,.  f dr 

0 


7'^  (/•, , Z|  |0,  ) - J"  7)|z(r)  I exp 

where  z(0)  -=  z^,  and  z(/-|  ) “Z|, 

TM^r„.z„\r„_^,z„_\  ) - J 7)|z(r)l  exp 
withz(r„_|)  —z„_|,and  z(r^)  ■=  z„;  and 


+ z (r)M  (r) 


(6.53a) 


/ dr 

'r  I 


+ zA/l 


(6.53b) 


R 

. 

2 

J"  D[z(r)  1 exp 

d^o  S dr 

'n 

1 

2 \dr 

+ zM 

(6.53c) 


with  z(r^)  “ z„,  and  z(/?)  — z.  These  expressions  differ  only  in  the  endpoint  conditions 
satisfied  by  the  paths.  We  now  focus  on  the  integration  over  M which  appears  in  Eqs.  (6.53). 
For  Eq.  (6.53a)  we  have 


J dr  z(r)M(r)  — — ^ kiZir^)  j"  drhir  — r, ) 


'o  1-1 


k\Z^r^  ) 
2k„ 


(6.54) 


since  the  r/s  are  ordered  and 


j*  drfiir  - r,  ) - y. 


Therefore 


T/^ir^,  Z|l0.  Zj)  - J"  £)(z(r)l  exp 


~2 


k 

. exp 

2ir/rj 

2 

'Fo(r,,zi|0,  Zj). 

0.  Further, 

+*«-|r.-l  )/2 

^ (r 

2) 


(6.55a) 


and 


giving 


^'M^r„.z„\r„_^.z„_^)  -e 

T^(R.z\r„.z„)  - z(r„,  z„), 

Tnf(R.z\0,z,)  ~ f d^^  </z„exp  j/  £ 

'*'o^r„,z„\r„_^.z„_^)  'F„(r,,z,|0,  2,). 
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Combining  Eqs.  (6.49),  (6.52),  and  (6.56)  and  using  Eq.  (6.51),  we  have 

oo  R 'n  '2 

'I'  («.  zlo,  z, ) - Y.  ( - Ik  f dr„  f r„_^  J dr^  f dz,  f dz„ 

n *0  0 0 0 

X (/(/■„,  z„)(/(r„_,.  z„_,  ) i/^r^,z^  ) 

X '•'„(«,  zl/-„,  z„)'»'„(r„,  zj z„_,  ) '^„(/-,,  Z||0.  z,),  (6.57) 

where  we  have  used  Eq.  (6.46)  to  do  the  integrals  over  the  /(,’s.  This  expression  is  the  standard 
Neumann  series  which  expresses  4'  as  a power  series  in  U.  If  Mr  is  approximated  by  keeping 
only  the  first  n + 1 terms,  one  is  using  nth-order  perturbation  theory.  To  zeroth  order 

4r(/?,  z|0,  zj  - 'P„(/?,  z|0.  z,).  (6.58) 

If  two  terms  are  kept,  one  is  making  the  Born  approximation 

R 

(/?,  z|0,  Zj ) ” 4'^  (/?,  z|0,  Zj ) - ik„  ( dri  ( dzi 

0 

X 4'„  (/?,  z|  r, , Z|  ) (y(/-| , Z|  )4'„  (r, , Z|  |0,  z^ ).  (6.59) 

The  second  term  in  Eq.  (6.59)  is  called  the  Born  term. 

In  the  preceeding  we  have  expressed  4^  as  a power  series  in  U.  It  is  also  possible  to  do 
perturbation  theory  around  a part  of  U.  Suppose 

U(r.z)  ~U„(r.z)  + U'ir.z).  (6.60) 

where  U'  is  in  some  sense  small  and  the  parabolic  equation  could  be  .solved  if  it  were  absent. 
We  might  then  expand  that  part  of  the  exponential  it.  Eq.  (6.1a)  which  depends  on  U'  and 
proceed  as  before.  We  would  obtain 

R 

4^  (/?.  z|0.  Zj ) “ 4'y  (/?,  z|0,  z^ ) - /At,  r dr,  r dZ| 

"1) 

X 4' (R,  z|  r, , Z)  ) (/'(r, , z,  ) 4' {y  (r, , Z]  jo,  Zj ) (6.61) 

where  'V  u satisfies  the  equations 

- 2/A(„0  'I'jy  {r,  z|r'.  z')  - jd/  - 2k^,U„{r.  z)j  4' (r,  z|  r',  z')  (6.62a) 

and 

'P  (j  (r,  z|r,  z')  “■6(z  — z').  (6.62b) 

If  Eqs.  (6.62)  cannot  be  solved  exactly,  Eq.  (6.61)  still  might  be  useful  if  'Pi,  is  approximated 
using  some  other  technique,  say  ray  acoustics. 

Before  going  on  to  the  Rytov  approximation,  we  mention  one  minor  point  which  may 
have  given  trouble  to  some  readers.  In  writing  Eq.  (6.49),  we  moved  the  integrals  over  the  r's 
outside  the  integral  over  the  path  z(r).  By  use  of  the  discrete  representation  for  the  path  in- 
tegral, it  is  obvious  this  interchange  is  permitted. 

6.4  The  Rytov  Approximation 

In  the  previous  subsection  we  developed  an  expansion  for  4^  in  powers  of  U.  In  certain 
situations  it  is  preferable  to  use  a perturbation  scheme  based  on  an  expansion  of  logr4'  in 
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powers  of  U.  This  scheme  is  called  by  Fradkin  [12]  "modified  perturbation  theory."  If  in  this 
expansion  one  keeps  only  the  first  nontrivial  term,  one  has  made  the  Rytov  approximation 
174-76).  The  Rytov  approximation  essentially  amounts  to  exponentiating  the  born  term.  In 
this  subsection  we  show  how  it  may  be  obtained  using  the  path-integral  formalism. 


We  wish  to  evaluate  the  path  integral  in  Eq.  (6.49),  which  we  have  called  without 
ordering  the  r,  variables.  We  now  have  three  techniques  for  doing  this.  We  can  use  the  veloci- 
ty representation  and  perform  a translation  in  v (r),  we  can  use  the  discrete  representation  and 
carry  out  an  analysis  smilar  to  that  which  led  to  the  evaluation  of  the  path  integral  in  Eq. 
(6.1  Sb),  or  we  can  use  ray  acoustics  directly.  Ray  acoustics  amounted  to  expanding  the  action 
(eikonal)  A(z(r))  through  second  order  about  the  lay  trajectory  z*.  If  A is  only  quadratic  in 
zir),  which  will  be  the  case  if  U(r,  z)  is  only  quadratic  in  z,  independent  of  its  dependence  or 
r,  then  ray  acoustics  gives  the  exact  solution.  The  integral  Tn^  is  the  solution  for  4'  with 
il(r,  z)  — -zM(r).  Therefore  ray  acoustics  directly  gives  the  value  for  T^- 


Although  it  is  instructive  to  evaluate  using  all  three  techniques,  we  shall  only  record 
its  evaluation  using  the  velocity  representation.  We  have 

R 

(r)  1 8 (z  - Zj  - ^ dr\(r)  1 

X exp|/)to  J"  dr  + M(r)z^  + M(r)  J ds  v (s)  |.  (6.63) 


We  next  let 


v(r)  + v'(r), 

K 


giving 


ik  R z ^ z ^ ^ 

Tm  - C’‘P  —2 I ~ ^ (6.6Sa) 


H 

^ drv'ir) 

exp  ik„  ( dr 

0 

We  now  use  the  relation 


R t K 

f drM(r)  J ds\'{s)  ” ” J*  drv'{r)M{r), 

*0  "o  *0 


where 

r 

Mir)  -^dsM(s). 

Equation  (6.66)  follows  from  the  constraint 


R 

^ drv'ir)  — 0. 
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Substituting  F.q.  (6.66)  into  7^  and  dropping  the  prime  on  v'(r)  gives 


H 

R 

f drvir) 

exp 

ik„  J dr 

0 

T) 

4"  V ^ (r)  — M(r)v  (r) 


Hguation  (6.68)  is  easily  evaluated  using  the  transformation 

V (r)  = wir)  + M(r). 

We  find 


Tm  = 


_kn_ 

2niR 


! 1/2 


exp 


ik..R 


jfdrM(r)  -±fdrM^(r) 


(6.68) 

(6.69) 

(6.70) 


Collecting  expressions  yields 


J D[z(r)  1 exp 


ik„  j dr 


+ z{r)M(r) 


= 'l'„  (R,  z|0,  z^.)  exp 


n 

/A„  J*  dr  Mir) 


X exp 


,k„R 

~ 2~ 


^ 1 " 


f dr  Mir)  - -|r  r dr  M^ir) 
o *0  7< 


(6.71) 


where  M is  defined  by  Eq.  (6.67).  We  have  thus  evaluated  another  important  path  integral. 


With  Mir)  defined  by  Eq.  (6.48)  we  find 

R 

r dr  Mir)  z,  + — - (z  - z^ ) 

7)  A 


I 


' - « 


O / -1 

r 


1 " 

r, 

= f LA- 
'^0  / « 1 

and 


~ j dr  M ir)  J dr  M i r) 
i)  i> 

R R A 

J dr  iMir))^  = f ds  Mis)  f ds'Mis')  f dr  f^ir  — s)Hir  — s') 
o *ri  n I)  A Ti 


(6.72a) 

(6.72b) 


I 


■f,,  v-i 


max  (r,.  r, ) 
' - — 


(6.72c) 


so  that 


^ 1 ^ 


ffdrMir)  -±JdrM^ir) 


1 


ik„R)^  : 


Y.k^r,iR  - r,) 


2 


ik„R)^  “ 


^ A,A^  |r/^  + 7?  max  (r,,  Zy ) — Rir^ 


(6.72d) 
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Substituting  into  Tm  (Eq.  6.71)  gives 


Tm  z|0,z  J exp  ' Z +-^  (z  -Zj)  " y ^ 


<>/ 


One  can  check  the  validity  of  Eq.  (6.73)  by  noting  that  if  the  r,’s  are  ordered,  so  that  min 
(r,,  r^)  = rj(i  > j),  Eq.  (6.73)  must  agree  with  Eq.  (6.56).  This  check  is  a straightforward  ex- 
ercise which  is  easily  carried  out  by  Fourier  transforming  with  respect  to  the  depth  variables  all 
the  'F’s  which  appear  in  Eq.  (6.56).  The  integrals  over  the  z,’s  then  yield  8 functions  which 
permit  all  but  one  of  the  integrations  over  the  Fourier  transform  variables  to  be  carried  out. 
The  last  integration  is  evaluated  using  Eq.  (4.21).  The  final  result  agrees  with  Eq.  (6.73). 

The  Rytov  approximation  follows  from  ignoring  in  the  cross  terms  involving 
(/  j).  Dropping  these  terms  and  substituting  Eq.  (6.73)  back  into  (6.49),  we  obtain 


S' (/?.  z|0,  z, ) = S',//},  z|0,  Zj)  Z dr„  { 

n -0  ”•  0 

X f/(r,,k,)  U(r„.k„) 


(Itt  ) (2ir  ) 


«P  -I  »,  + -2,)  -jj;'.  1 


-f.cs.rlo.z.)  i i 


X exp  ik  z,  + 


in  , / 

k:'' 


The  2n-dimensional  integrals  factor  into  n identical  two-dimensional  integrals,  which  permits 
the  infinite  series  to  be  summed  giving  the  simple  expression 


S' (7?,  z|0,  Zj ) •“  S*!,  (R,  z|0,  z ) exp  A", 


(6.75a) 


where 


K 

^ 'ko  S S 


U(r.  k) 


■I  . / \ ik^  . r 

t,p  ,k -—,1  -- 


(6.75b) 


Before  proceeding,  we  note  that  if  we  drop  the  k^  term  in  Eq.  (6.75b),  we  recover 
straight-line  geometric  optics: 


R 

- ik„  f dr  U r.  2, + -^  (z  - z,)  . 
0 K 
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Equation  (6.7Sb)  for  A' may  be  simplified  by  using 

U(r,k)  = J dz^e  Zj) 

and  performing  the  integration  over  k with  the  help  of  Eq,  (4.21): 

k.R 


(6.77) 


1/2 


X exp 


ik„ 


2/-,  R{R  - r,  ) 


2»r//|  {R  - /•,) 

[/?  (z,  - Z|  ) + r,  (z  - Zj ) 1 ^ 


(6.78) 


Since 


[/((Zj-Z|)  + Zi  (z  -Zj)l^  - R[R  - /-I  ) (Z|  -Zj)^ 

+ Z|/?(z  — Z|)^  — rj  (R  — rj  ) (z  —z^)^,  (6.79) 


we  have 


" VjR'r^tzJ  { ^ ')'„(/?.  zlr,.z.)(/(r,.z,)'|/„(r,.z,|0,zj.  (6.80) 

giving  the  Rytov  approximation  for  the  field  H': 


H'  (R,  z(0.  z, ) = 'I',,  (/?,  z|0,  z, ) exp 


- ik. 


, I { dn  f 

!'„(«.  z|0,z,)  *5  ' 

X (R.  z|  r, , z,  ) (/(r, , z,  )'R„  (r,  z,  jO.  z^) 


dz. 


(6.81) 


Comparing  Eq.  (6.81)  with  the  Born  approximation,  Eq.  (6.S9),  we  see  that  in  the  Rytov  ap- 
proximation the  field  has  the  form 


'R  = exp 


Born  term 


(6.82) 


The  Rytov  approximation  may  be  generalized  by  breaking  U up  according  to  Eq.  (6.60) 
and  treating  (/'as  a perturbation.  The  final  result  is 

[ modified  Born  term 


'R  “ 'R  (y  exp 


(6.83) 


where  'V  i,  is  given  by  Eqs.  (6.62)  and  the  modified  Born  term  is  the  second  term  on  the 
right-hand  side  of  Eq.  (6.61). 


Equation  (6.83)  is  an  important  result.  If  is  the  contribution  due  to  the  mean  sound 
speed  and  (/'  represents  the  effect  of  random  ocean  variability  on  the  sound  speed,  then  Eq. 
(6.83)  gives  an  adequate  description  of  the  propagation  process  over  most  of  the  region  in 
range  and  frequency  where  one  is  interested  in  doing  passive  listening.  Moreover  F'.q.  (6.83)  is 
amenable  to  the  direct  calculation  of  various  coherence  functions.  However,  one  would  prob- 
ably want  to  include  cross-range  dependence  by  using  the  three-dimensional  parabolic  equa- 
tion, and  the  numerical  integrations  which  would  be  necessary  are  not  trivial. 


The  work  in  this  section  points  out  what  seems  to  be  a general  characteristic  of  Feynman 
path  integrals.  High-frequency  approximations,  where  some  type  of  ray  picture  emerges,  are 
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natural  to  derive  using  path  integrals.  It  is  rather  awkward,  on  the  other  hand,  to  use  the  for- 
malism to  develop  perturbation  schemes. 


7.  BOUNDARY  CONDITIONS 

So  far  we  have  assumed  that  the  medium  is  infinite.  Now  we  come  to  the  crucial  problem 
of  incorporating  into  the  analysis  realistic  boundary  conditions.  Although  our  general  under- 
standing of  this  problem  is  rather  unsatisfactory  at  present,  the  path-integral  formalism  does 
give  some  indication  of  the  direction  future  research  might  take. 

The  goal  here  is  to  use  the  formalism  to  develop  algorithms  for  numerically  solving  the 
parabolic  equation.  With  the  inclusion  of  boundaries,  many  of  the  previously  derived  expres- 
sions do  not  apply,  and  it  is  almost  hopeless  to  obtain  even  an  approximate  solution  without 
detailed  computer  calculations.  We  are  always  assuming  the  medium  possesses  acoustic  charac- 
teristics which  are  range  dependent.  If  these  characteristics  arc  not  present  or  can  be  ignored, 
then  there  is  no  compelling  reason  for  introducing  the  parabolic  equation  in  the  first  place. 

The  boundaries  to  be  considered  are  the  surface  of  the  ocean  and  its  bottom.  The  surface 
of  the  ocean  does  not  create  any  particular  problems,  because  it  is  almost  always  assumed  in 
long-range-propagation  modeling  to  be  a free,  plane  surface  which  does  not  vary  with  time. 
When  the  ocean  bottom  is  considered,  a distinction  must  be  made  between  propagation 
through  water  with  a depth  excess  and  bottom-limited  or  shallow-water  propagation.  In  the 
case  of  a depth  excess,  the  exact  boundary  conditions  imposed  on  the  pressure  at  the  bottom 
interface  are  not  so  important  over  long  propagation  distances  because  acoustic  energy  which 
interacts  with  the  bottom  is  scattered  or  absorbed  and  can  for  the  most  part  be  ignored.  (This 
statement  requires  some  modification  if  the  receiving  hydrophone  or  array  is  bottom  mounted.) 
The  bottom  can  therefore  be  modeled  phenomenologically  as  being  perfectly  absorbing,  and 
discontinuities  in  density  and  sound  speed  can  be  ignored.  For  such  a situation  it  is  relatively 
easy  to  olve  the  parabolic  equation  using  the  split-step  Fourier  algorithm.  The  only  range- 
dependent  quantity  of  interest  is  the  sound  speed. 

For  cither  shallow-water  propagation  or  propagation  through  bottom-limited  water,  the  si- 
tuation is  far  more  complicated.  The  bottom  has  a pronounced  effect  on  the  received  signal  and 
cannot  be  modeled  as  casually  as  in  the  case  of  a depth  excess.  One  must  deal  realistically  with 
the  acoustic  parameters  of  the  bottom  material,  including  ab.sorption  and  discontinuities  in 
density  and  sound  speed.  Moreover,  variations  in  range  of  the  water  depth  and  of  the  bottom 
material  become  important,  as  well  as  the  range  variations  of  the  speed  of  sound  in  the  water. 
There  does  not  seem  to  be  an  algorithm  for  solving  the  parabolic  equation  in  this  situation.  All 
we  can  do  in  this  review  is  indicate  how  path-integral  techniques  may  be  used  to  develop  one. 


7.1.  Water  with  a Depth  Excess:  The  Split-Step  Fourier  Algorithm 


If  the  acoustic  energy  is  constrained  to  a channel  of  depth  L,  then  the  composition  law 
takes  the  form 

L 

(r,  z| /•',  z')  — ^ r/z, 'i' (r,  z|s,  z, )'!' (j,  z,|r'.  z').  (7.1) 
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II  is  easy  lo  argue  physically  for  the  limits  of  integration  in  Eq.  (7.1).  The  composition  law  is 
simply  the  mathematical  statement  that  the  propagation  process  can  be  viewed  as  occurring  in 
two  steps.  Sound  propagates  first  from  the  source  at  (r',  z')  to  an  intermediate  point  (s,  z^). 
The  function  'J' (s,  zjr',  z ) is  the  source  strength  of  a secondary  sound  field  which  originates 
at  the  intermediate  point  and  propagates  on  to  the  receiver  at  (r,  z).  The  field  at  (r,  z)  is  the 
coherent  sum  of  the  contributions  from  all  possible  intermediate  points  at  the  particular  range 
coordinate  s.  Consequently,  to  determine  the  limits  of  the  integration  in  the  composition  law, 
one  need  only  determine  those  points  which  could  serve  as  the  coordinates  of  an  acoustic 
source  (those  depth  positions  where  one  could  conceivably  place  a source).  For  the  channel 
being  considered,  a source  could  be  placed  anywhere  in  the  range  0 < z^  < L. 

If  we  take  the  intermediate  point  close  to  the  receiver,  we  have,  with  an  obvious  change 
of  notation, 

/. 

Mr  (/•  + Ar,  z|  €,  Zq  ) = f ifz'  Mr  (r  + Ar,  z|  r,  z'lMr  (r,  z'l  €,  Zq  ).  (7.2) 

*0 

Any  marching  algorithm  which  determines  the  field  at  r + Arfrom  the  field  at  rwill  necessari- 
ly involve  finding  some  expression  for  Mr  (/•  -I-  Ar,  z\r,  z')  and  then  evaluating  the  integral  over 
z'in  Eq.  (7.2). 


We  previously  concluded  that  for  Ar  sufficiently  small  ^ (r  + Ar,  z|r,  z')  is  approximately 
equal  to  the  solution  to  the  parabolic  equation  with  constant  index  of  refraction.  That  is. 


Mr  (r  -I-  Ar,  z\r,z') 


= + Ar,  z\r,  z')  exp 


2 


[n^  (r,  z) 


(7.3) 


where  Mr^  is  the  solution  to  the  parabolic  equation  with  n{r,  z)  = 1: 

— „(.r,  z\r',  z')  = 3?Mr„  (r,  z|  r',  z')  (7.4a) 

and 

Mr„(r',  zlr',  z')  = 6 (z  - z').  (7.4b) 

Both  fields,  Mr  and  Mr„,  satisfy  the  same  boundary  conditions  at  z = 0 and  z = L These 
conditions  are  in  turn  the  same  as  those  satisfied  by  the  acoustic  pressure.  To  derive  the  split- 
step  Fourier  algorithm,  one  assumes  that  L is  the  depth  of  the  water  column  plus  an  artificially 
introduced  bottom  layer  having  a depth  1/4  to  1/3  the  depth  of  the  water  and  that  the  field 
vanishes  at  z — 0 and  at  z — L. 


and 


'V„{r.O\r',z')  -0 
L\r'.2')  -0. 


{7.5a) 

(7.5b) 


Equations  (7.4)  and  (7.5)  are  easy  to  solve.  One  expands  the  field  Mr^  in  terms  of  a complete 
set  of  eigenfunctions  which  satisfy  Eqs.  (7.5): 


'VJr.zlr'.z')  - 1 r'.  z') 


2 

1/2 

knz 

L 

sin 

L 

k -1 

Substituting  Eq.  (7.6)  into  Eq.  (7.4a)  and  using  the  orthonormality  property 

L 

sin  I sin 


-r  j dz  sir 
L i 


k'nz 


^kk‘ 


(7.6) 


(7.7) 
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we  get 


or 


2ikgb,a^  (r,  r',  z') 


kiT 


a^(r,r'.z')  - a*  (r',  r',  z')  exp 
Since  we  have  the  completeness  property 


a*  (r,  r',  z'), 


iir  — r')  I kn 
L 


2k„ 


2 

^ . knz 

sin 

■ 

knz' 

L 

L S'"  I 
* -1  ^ 

L 

- 8 (z  - z'). 


Eq.  (7.4b)  will  be  satisfied  if 


a*(r',r',  z')  - — 


1/2 


sm 


nz'k 


(7.8a) 


(7.8b) 


(7.9) 


(7.10) 


independent  of  r'.  Combining  Eqs.  (7.3),  (7.6),  (7.8b),  and  (7.10),  we  find 

[ i^rk„  . 1 

(r  + ^r,  z\r,  z')  - exp  | — - — ln^(r,z)  - 1 )[ 


X S'" 


k -1 


kiTZ 


sm 


knz' 


exp 


2k 


i^r  \ kn 
L 


(7.11) 


Substituting  Eq.  (7.11)  into  Eq.  (7.2)  and  interchanging  the  integration  over  z' with  the  sum 
over  k,  we  get 

'I' (r  + Ar,  z|e,  Zo)  -expj — - — (n^(r,  z)  “'1} 


1^1^  . knz 


exp 


/Ar 

kiT 

2ko 

L 

’'I*' 


sin  S' (r,  z1«,  Zq  ). 


(7.12) 


When  one  assumes  the  field  vanishes  at  z “ L,  one  introduces  undesirable  reflections  off 
the  artificial  bottom.  They  are  eliminated  by  adopting  the  point  of  view  that,  in  water  with  a 
depth  excess,  the  bottom  acts  as  a perfect  absorber  and  any  acoustic  energy  which  strikes  the 
bottom  never  gets  back  up  into  the  water  column.  The  simplest  way  of  mathematically  realiz- 
ing this  notion  is  to  add  an  imaginary  term  to  the  index  of  refraction: 

n^ir.z)  ^ n^(r.z)  + i^iz).  (7.13) 


where  ^ is  small  until  z is  close  to  L.  It  is  typically  taken  to  be  a Gaussian; 

, I I . 1 21 


((z) 


exp 


Lt 


(7.14) 
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Therefore 


[nUr.z)  - 11 


exp  — Ar/>  exp  — 


\ L — z 
I L*r 


ln^(r.  z)  - 1 1 


Since  the  inclusion  of  ( is  phenomenological,  the  parameters  p and  <t  can  be  determined  only 
through  experience  gained  by  examining  the  numerical  results  for  various  test  cases. 

The  next  step  is  to  approximate  the  integral  over  z'by  a discrete  sum.  It  is  assumed  that 
over  a depth  increment  Az  the  field  may  be  taken  to  be  constant.  Let  A/ be  the  number  of 
increments:  N^z  —L.  Then 


r ^ 

f dz'~-  Az  X 

0 m — I 

and  Eq.  (7.12)  takes  the  form 

'P  (r  + Ar,7Az|e,  Zq  ) “ exp  I —Arp  exp  ~ — 


(ArA„  , 

exp  — [n^ir.j^z)  -1) 


„ 2 ^ . \kTTj\  -/Ar  \ kn 

X T sin  'i>  (r.  m£kz\t,  Zf.).  (7.16) 

In  writing  Eq.  (7.16)  we  have  used  Eq.  (7.15)  and  have  truncated  the  sum  over  k.  This  truncat- 
ed is  consistent  with  the  assumption  that  the  wiggles  in  'I'  of  spatial  size  less  than  AZ  are  unim- 
portani.  It  is  analogous  to  the  usual  Fourier  sampling  theorem.  Although  the  sums  in  Eq. 
(7.16)  go  from  1 to  A,  we  could  sum  from  0 to  A — 1 or  from  I to  A — 1 without  changing 
the  results. 

So  far  we  have  assumed  one  is  inlereited  in  calculating  the  Green’s  function  in  which 
case  'Ke,  wAzle.  zq)  - /Az.  where  zo  - w„Az.  If  one  is  interested  in  calculating  the  solu- 
tion i/((r,  z)  to  the  parabolic  equation  corresponding  to  some  initial  distribution  h, 

i/<  (f,  z)  — h(z), 

then  Eq.  (7.16)  is  replaced  by 

i/i  (r  -E  Ar.yAz)  —exp  — Arp  exp  — exp  | (n^Ir.yAz)  ~ ' l| 


2 

N 

G t n 

knj 

— (Ar  kn 

N 

^ Sin 
k -1 

2k„  N:kz 

where 


X ^ sin  il)(r,mAz), 

m-l  ^ 


i/i(«,mAz)  — fi(mSz). 
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Let  US  examine  the  structure  of  Eq.  (7.17)  by  first  simplifying  the  notation.  We  define, 
for  j - 0 N, 


+ Ar.yAz)  — ^j{r  + Ar). 
^(r.jAz) 

h(JAz)  - hj, 


exp 


exp 


— Arp  exp  - 
iArk_ 


Ntr 


-«>• 


In^fr.yAz)  - 111  -&j(r). 


and 


Therefore 


2 

— /Ar  jir 

A 

2k„  NAz 

-^r 


4>j{r  + Ar)  ~Qj&j(r)  X ^ L 

* -0  m -0 


</<„(r). 


(7.19a) 

(7.19b) 

(7.19c) 

(7.19d) 

(7.1 9e) 

(7.190 

(7.20) 


This  equation  defines  the  split-step  Fourier  algorithm.  It  is  easy  to  see  from  Eq.  (7.20)  how  the 
programming  should  go.  To  calculate  the  field  at  r -I-  Ar,  one  does  a sine  transform  on  the  field 
at  r.  The  result  is  then  multiplied  by  9^,  and  a second  sine  transform  is  performed.  The  output 
of  this  second  transform  is  multiplied  by  3j  times  6y(r),  giving  ^j(r  + Ar).  The  pressure 
is  obtained  by  multiplying  (r  + Ar)  by  some  constant  times  r'^^exp  ik„(r  -I-  Ar)  (see  Eq. 
(3.4)).  After  the  value  of  the  pressure  is  transferred  to  a storage  location,  the  process  is  repeat- 
ed until  the  maximum  range  of  interest  has  been  reached.  One  initializes  by  taking 

-h„.  (7.21) 

The  matrices  3^  and  7,  are  calculated  once  and  stored.  The  matrix  (r)  depends  on  the  varia- 
bility in  range  and  depth  of  the  sound  speed  and  is  therefore  recomputed  at  each  range  step. 


If  one  does  not  have  a subroutine  for  doing  a fast-Fourier  sine  transform,  the  usual  ex- 
ponential FFT  may  be  used  by  doubling  the  space.  The  expression 

N-\  («.!.  1 JA-I 

2/  X Sin  fk  - 1.  (7.22) 

A .0  I ^ I * ”0 

is  used,  where is  defined  by  the  relations 
/o'*"  -/o  -0, 

*//v  -0. 

/f"  -/*.  * - 1 N - 1. 

and 

fdbie  . k ~N  + 1 2N  - 1.  (7.23) 
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Let  US  now  try  to  isolate  the  various  assumptions  which  went  into  the  derivation  of  the 
split-step  Fourier  algorithm: 

• The  parabolic  approximation  is  valid.  This  assumption  will  be  discussed  in  Section  9. 

• The  density  is  constant.  This  is  a valid  assumption  provided  propagation  through  the 
bottom  is  unimportant. 

• The  step  size  Arcan  be  chosen  small  enough  so  that  Eq.  (7.3)  is  valid.  As  long  as  the 
sound  speed  is  a smooth  function  of  depth  and  range,  a Arcan  always  be  found.  If  the  sound 
speed  possesses  discontinuities,  then  Eq.  (7.3)  is  questionable,  regardless  of  the  size  of  Ar  (and 
Az).  The  reason  is  rather  easy  to  see  physically.  Equation  (7.3)  means  that  the  ray  paths  follow 
straight  lines  from  r to  r + Ar.  If  the  sound  speed  is  discontinuous,  then,  in  addition  to  (or  in- 
stead of)  the  straight  line  path,  ray  paths  which  are  reflected  and  refracted  by  the  discontinuity 
will  contribute  to  the  field  at  r Ar.  These  reflected  and  refracted  paths  are  not  included  in 
Eq.  (7.3).  The  problem  obviously  persists  as  Az  — > 0. 

• In  water  with  a depth  excess,  acoustic  energy  which  interacts  with  the  bottom  is  ab- 
sorbed or  scattered  and  does  not  contribute  to  the  received  signal.  There  is  some  experimental 
evidence  in  favor  of  this  assumption  [77].  For  a source  close  to  the  surface  it  will  introduce 
some  bias  in  the  calculated  pressure  in  the  regions  between  convergence  zones  but  probably 
not  in  the  neighborhood  of  a convergence  zone,  where  RR  and  RSR  paths  are  so  much  more 
important  than  BRSR  paths. 

• The  preceding  assumption  may  be  mathematically  modeled  according  to  Eqs.  (7.5b), 
(7.13),  and  (7.14).  One  technique  is  probably  as  good  as  another.  There  is  a situation,  however, 
where  care  should  be  taken.  If  one  is  numerically  investigating  the  validity  of  the  parabolic  ap- 
proximation (or  some  improved  approximation)  by  comparing  the  output  of  a PE  program  witli 
that  of  a normal-mode  program,  one  may  be  considering  differences  of  just  a few  decibels.  In 
such  a situation  one  must  be  sure  the  absorbing  bottom  is  modeled  identically  in  both  pro- 
grams. Otherwise  the  results  may  not  really  say  much  about  the  validity  of  the  parabolic  ap- 
proximation, particularly  in  the  regions  between  convergence  zones.  As  an  example,  one  would 
a priori  question  the  conclusions  of  a comparison  of  the  field  calculated  using  Eq.  (7.20)  with 
that  calculated  using  a normal-mode  program  if  the  absorbing  bottom  is  modeled  in  the 
normal-mode  program  by  cutting  off  the  modal  sum  when  the  phase  speed  exceeds  the  sound 
speed  at  the  bottom.  It  would  be  preferable  to  compare  the  output  of  a program  based  on  Eq. 
(7.20)  with  a modal  solution  where  the  normal-mode  eigenfunctions  Z„(z)  are  calculated  using 


Z„(z)  k^ZJz)  (7.24) 

with 

Z„(0)  ’-Z„(L)  -0.  (7.25) 

• A reasonable  value  for  Az  may  be  found.  The  choice  of  the  number  of  depth  points 
N “ L/Az  is  intimately  related  to  several  of  the  other  assumptions  in  a way  which  we  would 
now  like  to  discuss.  In  Eq.  (7.2)  we  assumed  Aris  small  enough  so  that  'F  (r  -F  Ar.  z|r,  z ) acts 


■f  kl  l«^(z)  -F  li  (z)l 
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almost  like  a 6-funclion  in  z — z'.  With  this  assumption  we  were  able  to  take  the  index  of  re- 
fraction to  be  a constant  and  hence  find  an  analytic  expression  for  M' (r  -I-  Ar,  z|r.  z').  The 
crucial  point  is  how  'I' (r  -I-  Ar,  z|r,  z')  imitates  the  behavior  of  a 8-function.  If  \z  — z'|  is 
large,  'I'  will  be  a rapidly  oscillating  function  of  z'.  Hence,  when  one  integrates  over  z',  cancel- 
lations will  occur  for  large  |z  — z'\  which  effectively  restrict  the  range  of  integration  to  a small 
region  about  z = z.  If  one  does  not  sample  often  enough  in  depth,  the  oscillations  are  washed 
out,  and  these  essential  cancellations  do  not  occur.  The  result  is  an  erroneous  value  for  the 
field  at  /■  -I-  Ar. 


On  the  other  hand,  one  condition  for  the  validity  of  the  parabolic  approximation  is  that 
the  slopes  of  the  Feynman  paths  should  be  small  (5):  | (z  — z')/Arl  1.  That  is,  the  contri- 
bution to  the  field  from  all  paths  with  | (z  — z')/Ar|  > 1 should  be  negligible.  This  small- 
slope  requirement  is  really  not  incorporated  into  the  numerical  algorithm.  If  the  algorithm 
could  be  modified  so  that  one  keeps  only  those  contributions  to  the  field  which  correspond  to 
1 (z  — z')/Ar|  <5C  1,  then  the  end  result  would  be  the  same  (if  the  parabolic  approximation  is 
valid).  More  importantly,  one  would  not  have  to  sample  as  often  in  depth,  because  the  oscilla- 
tions which  occur  for  |z  - z | large  would  no  longer  be  as  important.  Consequently,  one  would 
conceivably  save  storage  space  and  running  time. 


7.2.  Bottom-Limited  Water  and  Shallow  Water:  An  Unsolved  Problem 

Perhaps  the  most  important  unsolved  problem  involving  the  application  of  the  parabolic 
equation  to  sound  propagation  is  the  development  of  an  algorithm  for  obtaining  its  solution  in 
a bottom-limited  or  shallow-water  situation.  We  will  not  solve  this  problem  here  but  rather  in- 
dicate how  it  might  be  approached. 


For  concreteness  we  will  consider  a specific  model  for  the  acoustic  medium.  We  will  as- 
sume it  consists  of  two  fluids  in  contact  at  the  bottom  interface.  The  depth  of  the  water  is 
B(r),  and  the  thickness  of  the  bottom  fluid  is  L — B(r),  with  the  acoustic  pressure  being 
essentially  zero  for  depths  greater  than  L The  density  of  the  water  is  p^,  and  the  density  of 
the  bottom  fluid  is  Pf,.  Both  p^  and  Pf,  are  constant.  For  the  index  of  refraction  we  have 

n^ir.z)  z),  0 < z < Bir), 

“ (r,  z)  -f  /T)^  (r,  z),  Bir)  < z < L.  (7.26) 

where  represents  the  effect  of  bottom  absorption. 


If  the  water  depth  and  index  of  refraction  are  independent  of  range,  the  pressure  is  given 
by  a normal-mode  sum 

“TTTh' 

H p ^ 

where 

piz)  — P^,  0 < z < B, 


“ Pi,,  B < z < L, 

and  the  normal-mode  eigenfunctions  satisfy  the  differential  equation 


dz^ 


+ kln^(z) 


ZJz)  -kjlZ„(z). 


(7.28) 

(7.29) 
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the  boundary  and  continuity  conditions 

Z„(0)  ~ Z„(L)  -0,  (7.30a) 

Z„  (z)  continuous,  (7.30b) 

and 


1 d r,  I ^ 

— T": — — Z„  (z)  continuous, 
lt(z)  dz 


and  the  orthonormality  condition 


'■  // 


(7.30c) 


(7.31) 


These  functions  are  not  normalized  in  the  conventional  way  because  we  are  assuming  the 
pressure,  rather  than  the  particle  velocity  potential,  satisfies  the  Helmholtz  equation. 


In  Section  9 we  will  show  the  parabolic  approximation  is  equivalent  to  replacing  Eq. 
(7.27)  by 

' ,k  R 

e ” '\'{R.z\0.zJ  (7.32) 


p(x) 


1 

An 


2tti 


k„R 


with 


'ViR.  z|0,  z,) 


!>  ) 


£Z„(2)Z„(z,)  exp 


— (*2  - k^) 
2k  ' 


(7.33) 


By  use  of  the  orthonormality  condition  (Eq.  (7.31)),  it  is  easy  to  show  that  9'  obeys  the  compo- 
sition law: 


i. 

M' (/•  -I-  Ar,  z|0,  z, ) = f rfz' 'P  (r  -I-  Ar,  z|  r,  z')'P  (r,  z'|0,  z.  ), 


(7.34) 


where 


Mr  (r  -I-  Ar,  z]  r,  z') 


1 

It  (z') 


5^  Z„(z)Z„(z')  exp 


(k^ 


(7.35) 


With  Eqs.  (7.34)  and  (7.35)  we  have  a marching  algorithm  for  solving  the  parabolic  equa- 
tion for  the  range-independent  situation.  Of  course  one  would  not  use  these  equations  in  prac- 
tice since  the  pressure  is  readily  obtained  from  Eq.  (7.27).  V^e  have  derived  them  to  fix  the 
notation  and  because  their  structure  suggests  how  one  might  generalize  to  the  range-dependent 
problem. 


As  we  have  mentioned,  the  problem  of  constructing  an  algorithm  in  the  general  case 
reduces  to  the  problem  of  determining  an  expression  for  'P  (z  -I-  Az,  z|z,  z ).  If  one  has  such 
an  expression,  the  composition  law  can  be  used  to  construct  the  field  at  z -f  Azin  terms  of  the 
field  at  z.  If  the  parabolic  approximation  is  to  be  at  all  valid,  a value  of  Azcan  be  found  such 
that  the  variation  with  range  of  the  medium’s  acoustic  parameters  can  be  ignored  within  the 
region  z to  z -t-  Az.  With  this  approximation  'P  (z  -f  Az,  z]  z,  z')  has  the  form  of  Eq.  (7.35)  but 
with  the  usual  normal  modes  replaced  by  range-dependent  normal  modes  Z„(z,  z)  [78].  These 
range-dependent  eigenfunctions  and  their  corresponding  eigenvalues  are  obtained  by  solving 
an  eigenvalue  problem  at  the  particular  range  z. 
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Specifically  we  have 

H'  (r  + Ar,  z\  r.  z')  - (r,  z)Z^  (r,  z')  exp  [k^  (r)  - k},  ]|,  (7,36) 

where 

Z„(r.z)  - Z*(r,  2).  0 < z < B{r). 


- Z^(r.z).  B(r)  < z < L. 

(7.37) 

with 

^ +kf,nl(r.z)  Z;(r.z)  ~ kl  (r)Z”  (r.  z) 
dz^ 

(7.38a) 

for  0 < 

z < B(r)  and 

+ /c,^  (Wp  (r,  z)  +/T)^(r,  z))  Z^  (r,  z)  k^  (r)Z^  (r,  z) 

\Bz^ 

(7.38b) 

for  0(r) 

< z < L.  The  boundary  and  continuity  conditions  are 

Z,‘^(r.O)  ^Zl^U.L)  -0, 

(7.39a) 

Z;(r.  fl(r))  ~ Z!;(r.  Bit)). 

(7.39b) 

and 

(7.39c) 

If  we  further  define 

Mr  (r,  z|o.  z,)  z),  0 < z < Bir). 

- '1'  *(2,  z),  Bir)  < z < L, 

(7.40) 

then  the  pressure  at  a point  (R,  z)  in  the  water  column  is 

( I 1/2 

e'*”" z). 

where  'V  is  obtained  by  using  the  algorithm 

'f'  "'(/■  + Az,  z)  ~e  £ z„”'(r,  z)e'^''‘" 

fl 

1 1 y 

X j-^  I rfz'Z„*'(r,z')'«'"'(r,  z') 

+ — f rfz'Z.Mr,  z')'l'*(r.  z')|  (7.41a) 

Pa  rfir)  ) 
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and 


'P 


Ur  + Ar.  z)  -=  Z*  (r.  z)e'^"’' 


B(r) 


l>w  *& 
1 '• 


J </z'Z„"'(r,  z')M'‘"(r.  2 ) 


f>h  B(r) 


I 


dz'  Z^(r.  z')'l'*(r,  z') 


(7.41b) 


Since  this  approach  has  never  been  tested,  it  is  not  appropriate  to  discuss  the  algorithm, 
Eqs,  (7.41),  in  detail.  We  would  like  to  make  two  points  however.  First,  we  have  ignored  any 
contribution  due  to  shear  waves.  It  is  possible  in  principle  to  generalize  the  whole  formalism  to 
include  shear  waves  in  the  same  way  quantum  mechanics  is  generalized  to  include  particle  po- 
larization (spin). 


Second,  range-dependent  normal  modes  have  been  used  (78)  to  solve  the  Helmholtz 
equation.  The  main  problem  with  their  use  is  that  the  coupled  differential  equations  for  the 
modal  coefficients  are  extremely  difficult  to  solve.  Equations  (7.41)  really  amount  to  solving 
these  equations  in  the  parabolic  approximation.  However  Eqs.  (7.41)  repre,serit  an  unsophisti- 
cated, brute-force  technique  for  solving  the  parabolic  equation.  A complete  .set  of  normal 
modes  mist  be  calculated  for  each  range  step.  One  imagines  there  is  a simpler,  more  sophisti- 
cated technique  than  one  based  on  range-dependent  normal  modes. 


8.  PATH  INTEGRATION  AND  THE  HELMHOLTZ  EQUATION 

In  this  section  we  will  develop  path-integral  representations  for  the  solution  to  the 
Helmholtz  equation  and  consider  a few  simple  applications.  We  shall  assume  the  medium  is 
unbounded.  It  will  be  apparent  that  many  of  our  results  are  valid  even  for  a bounded  medium. 

8.1.  Some  Basic  Representations 

We  begin  with  the  Helmholtz  equation,  Eq.  (3.2): 

(V2  -t-  k2o2(x)lp(x) 6'J»(x  -X,).  (8.1) 

where  — w/c„,  n(x)  — f„/f(x),  and  we  have  suppressed  any  time  dependence.  With  an 
integration  by  parts,  one  can  show 

p(%)  e"  ■^''’^^<I>(t,x|0.  X,).  (8.2) 

where  the  propagator  <t>  satisfies  the  equations 

-2/*„d,«D(T,x|T',  x')  - (V^  -E  *„M«^(x)  - 1 1|«I>(t,  x|t',  x ) (8.3a) 

and 

<I>(t',  x|t',x')  -8‘^’(x  - x'),forallT'.  (8.3b) 

As  before,  it  will  be  convenient  to  define 

U(%) y |/»^(x)  - II.  (8.4) 
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Formally,  at  least,  we  have  the  composition  law 

‘h  (t,  x|0.  ) = J*  r/x ' <t>  (t,  xj (r,  x')<l>  ((r,  x'lO,  X, ),  (8.5) 

where  (t  is  specified  only  by  the  condition  t > o-  > 0.  If  we  iterate  Eq.  (8.5),  we  obtain 

r ^ 

<l>  (t,  x|0,  Xj  ) =J  rfxj  ..  rfXyv-l  n -I  >*/-)  )• 


with  (fg  “ 0,  Xg  = Xj.  ‘T  — T,  and  Xyy  = X.  We  will  choose  increments  of  equal  length, 

<T , — <r  I = Ait,  (8.7) 

for  all  I.  We  now  assume  N h large  enough  (Ait  is  small  enough)  so  that  U(x)  can  be  con- 
sidered constant  when  calculating  each  of  the  component  propagators  in  Eq.  (8.6).  We  there- 
fore have 


4>(o-,,  x,|<r,_,,  x,_|  e (<r,.  x,|rr , , x,  _,  ), 

where  satisfies  the  equations 

- 2/7c„d„  d>„(o,,  x,la,_,,  x,_,  ) = S7}<1>Jit,.  x,1ci-,_,  , x,_,  ) 
and 

4>„(»r,_t.X,trr,_i.x,_,  ) =6'^>(x,  -X,_,). 

The  solution  to  Eqs.  (8.9)  is 

'F„(»i-,,x,lo-,_,,x,_|  ) 

Gathering  expressions  gives 


ko 

3/2 

ik„AiT 

X,  - x,_, 

2 

Ini  Ait 

exp 

2 

Atr 

(8.8) 

(8.9a) 

(8.9b) 

(8.10) 


d*  (t,  x|0,  X, ) = 


Ini  Ait 


3/V/2 


X exp 


ik„AiT  j; 


X,  -x,_, 


Ait 


- i/(x,) 


(8.11) 


In  the  continuum  limit,  yv—  <»,  with  NAit  = t fixed,  Eq.  (8.1 1)  becomes  exact.  We  have 

IXg,  X|  , • ■ • , Xyyl  — X (fr  ).  0 ^ (r  < T, 


*1  - *0  *2  - *1  *\  - */v-i 


Ait 


Ait 


Ait 


d\  (it  ) 

diT 


Ait  ^ f diT, 

I -1  0 


and 


|3A/2 


IttiAit 


f d\^  . <fX;v_l  — J"  Dlx((r)l, 


(8.12a) 

(8.12b) 

(8.12c) 

(8.l2d) 
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giving 

<I>  (t.  x|0,  Xj ) = J*  D(x  (cr  ) I exp  (A„  J"  iAr  y — -U<x(ir))  . (8.13) 

According  to  Eq,  (8,12a)  the  paths  x(rr)  in  Eq.  (8,13)  satisfy  the  endpoint  conditions 

x(0)  “ Xj.  x(t)  = X.  (8.14) 

Substituting  Eq.  (8.13)  into  Eq.  (8.2)  gives 

p(x)  = ^ ? r/T  J D[x(/-)1  exp  (/c„  j”  rf<r  y -U(xUr))  . 

2k„  i •'  0 M " 

(8.15) 

Because  the  measure  depends  on  t,  the  integral  over  t cannot  be  interchanged  with  the  integral 
over  the  paths.  Using  Eq.  (8.4)  we  can  rewrite  Eq.  (8.15)  in  the  form 

2 

p(x)  = — 1 dr  J DlxUr)]  exp  ^ J dtr  + n^ixUr))  . (8.16) 

2k^  i 2 “5  d(r 

where  again  the  integration  is  over  all  continuous  paths  satisfying  the  endpoint  conditions,  Eq. 
(8.14).  This  is  our  first  path-integral  representation  for  the  solution  to  the  Helmholtz  equation. 

Just  as  with  the  parabolic  equation,  it  is  possible  to  introduce  a change  of  variable 

^ , (8,7) 

and  construct  the  velocity  representation.  Taking 
(.  , \3/V/2 

\-^\  J dxy  ...  J DixUr)].  (8.18) 

we  have 

T 

<t»(T,  x|0,  Xj)  “ J £>[v(fr  ))8 X - Xj  -^dirvitr) 

X exp  |/ko  J"  </rr  -yV^fcr)  — (/|xj  + j"  dcr  ' X ((T  ) (8.19) 


We  therefore  have  our  second  representation  for  the  pressure  field; 
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« r j 

p(x)  = — ^ J*  </t  r Z)(v  (rr  ) J6  X — X^  — f d(r\  (tr  ) 
2^0  0 *0 


X expj-^  rfrr  |v  ^ ((r  ) + + J citr’yUr  ) (8.21) 

Apparenily  this  particular  represenlalion  was  first  obtained  by  Fradkin  112). 

We  shall  construct  a third  representation  by  introducing  a new  computational  technique. 
Basically  we  want  to  redefine  the  path  x(<r),0  < <r  ^ t,  according  to  the  transformation 

w(r)  =-v/— x-^r  (8.22) 

and  integrate  over  w rather  than  x.  Here  R'  will  be  a free  parameter.  When  we  consider  appli- 
cations, we  will  let  /('be  the  range.  The  advantage  of  this  tran.sformation  is  that  the  resulting 
path  integral  will  have  a form  similar  to  that  of  the  path  integral  for  the  parabolic  equation. 

The  path  in  F.q.  (8.22)  is  defined  on  the  interval  0 < r < /(  ' and  satisfies  the  endpoint 
conditions 

w(0)  = y*(R')  = (8.23) 

For  the  exponent  in  Eq.  (8.13)  we  have 


i 2 ^ -(/(x(<r))  - S dr  \ 


-L  ^ 

2 dr 


(r)  . (8.24) 


For  the  measure  we  have 


( ^ 3m 


~ 2^1 

k R'  1 


,,  \3m 


R' 


where 


/ /;(x(.r)]  - f D(w(r)l, 


. 3 m 

l^r 

■ f 
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The  path  integral  in  Eq.  (8.13)  becomes 


•I*  (t,  x|0,  ) = 


3/2 


1’t/R  *|0.  *»). 


(8.28) 


where 


<!>,//(  (/(',  xjO,  *j ) — J"  DIw  (r)  1 exp 


R 

2 { /—  1 

ikf,  f dr 

X 

(M 

- -v/—  w(r) 

0 

2 

dr 

R'  1 V R'  \ 

(8.29) 


with  w(0)  •=  JTR'It)  X,  and  w(/?')  “ yfJWTr')  x.  Before  substituting  into  the  integral  over 
T which  determines  the  pressure,  we  let  fi  “ t/R'.  Hence 


p(x)  J (/?',  x|0,  xj. 


(8.30) 


2^„  i p- 

(It  is  implicitly  assumed  that  /(„  has  a small,  positive  imaginery  part.)  This  is  our  third 
representation  and  the  last  one  we  will  derive. 

8.2.  A Homogeneous  Medium 

If  the  medium  is  homogeneous  so  that  nix)  is  constant  ( = «„),  the  complicated  path- 
integral  representations  of  the  previous  subsection  must  give 

/*„  n„  I » — I 


pix) 


4w| X — xj 


(8.31) 


Let  us  verify  that  this  is  indeed  the  case  for  the  representation  Eqs.  (8.29)  and  (8.30).  It  is  not 
difficult  to  show  the  path  integral  in  Eq.  (8.29)  is  normalized  so  that 


/ Dlw(r)l 


exp 


ik 


n 


dm 


dr 


l<o 

3/2 

exp 

( 

m(R’)  - 

wJO) 

IttiR' 

1 2 

R' 

K 

3/2 

exp 

,k„R' 

* ~ P 

2mR‘ 

2P 

R'  \ 

(8.32) 


If  the  index  of  refraction  is  independent  of  x,  then 


R r 

jf  -/Jl/(^/i8w (/■))! R'PU  ~ ^ R'pin},  -1).  (8.33) 


We  therefore  have 

x|0.  X,)  - 


*0 

3/2 

,k,R' 

litiR' 

exp 

2 

X — Xj 

2 

R' 

(8.34) 
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Substiluling  this  into  Hq.  (8.30)  yields 
p{\)  = 


J 

ko 

\ /2  OO  i ■ 1 fi  * 

r dfi  'k„l< 

+ j 

X — 

2 

47T 

2ttiR' 

1 •»/'>  exp 

n 2 

By  letting  t = — *J,  this  equation  becomes 

p(x)  = 


1 

l/(  ” , 

r dr 

/A'„/I„|x  - x j 

r -eI 

47r 

27r/|  X — Xjl 

J 3/2 

0 T 

2 

T T 

T 

(8.35) 


(8.36) 


In  Appendix  A a one-parameter  integral  representation  for  the  free-field  Green’s  function  in  n 
dimensions  is  constructed.  By  comparing  Eq.  (8.36)  with  Eq.  (A  13),  we  obtain  the  desired 
result,  Eq.  (8.31).  Notice  that  /('canceled  out,  as  it  should,  since  Eq.  (8.30)  is  independent  of 

R'. 


8.3.  The  Parabolic  Proflle 

No  general  discussion  of  Feynman  path  integrals  would  be  complete  without  mention  of 
the  parabolic  profile  model.  In  this  subsection  we  shall  examine  this  model  and  point  out  the 
relationship  between  the  path-integral  solution  and  the  usual  normal-mode  solution. 

We  first  observe  that  if  is  a function  of  just  depth,  U(x)  = U(z),  the  path  integral  in 
Eq.  (8.29)  (with  fi  = t/R  ) may  be  written  as  the  product  of  three  path  integrals: 

(8.37) 


(/?',  xjO,  Xj)  - 


where 


and 


in  which 


<hv 


'h., 


f D[wy(r)]exp  f 


2 

ik., 


/* 

dWy 

2 

I 

- 

dr 

H 

dw., 

2 

f 

y 

0 

dr 

f Dlw^(r)}  exp  J ~ 2ik„fl  f drlKyf^w,) 


dw^ 


dr 


and 


(/?')  — fi 



7 

wJR')  -fi-'^l'^y. 

wJR')  -fi~'f^z. 

The  integrals  are  normalized  so  that 
i L \m 


27r(Arl 


J dw^^  dw,^^  ^ J Dlw Jr)]. 


(8.38a) 

(8.38b) 

(8.38c) 


(8.39) 
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with  similar  expressions  for  f D[w  (r)]  and  f Dltv,  (r)).  It  should  be  obvious  by  now  that 

I 1 1/2  I , 0. 1 ..  .121 


■*'  |2ir/7?'j 

Substituting  these  into  Eq.  (8.30)  gives 


-L  f-^q.exp  ^ + 1 . (8. 

477  Jo  ^3/2  2 ^ p R'  PR' 


X — Xj 

~2P  'F~ 


'*0^'  > ->1 


2p  R' 


1 2 2 
1 X - xj  . 1 y - 3', 


(8.40a) 


(8.40b) 


Although  we  could  keep  R'a  free  parameter  and  demonstrate  that  it  cancels  out  in  the  end,  it 
simplifies  the  notation  to  set  /?'-/?  ■=  l(x  -Xj)^  + (y  -yj)^]'^l  The  exponent  in  Eq. 
(8.41)  then  becomes 


ikgR  ] 

P + ~ 

2 p 


We  let 


= 1 - a2(z  - Zw)^ 


where  Zw  is  the  depth  of  the  SOFAR  axis.  Hence 


PU{y/pw^{r))  •=  y /3  ^ [Wj  (r)  -p 


Therefore 


•Pj  “ / ,1" -P^a^\^Ar)  - /3 '^^z^l^jj.  (8.44) 

We  know  that  ray  acoustics  gives  the  exact  value  for  the  path  integral  if  Uis  at  most  qua- 
dratic in  the  path  variable.  We  will  use  this  fact  to  evaluate  <IJj.  Since  <1Jj  has  exactly  the  same 
form  as  the  path  integral  we  considered  in  Section  3,  we  write 


|27r/dH'/(/?)/dp(0)  I 
where  p(0)  “ </w^*  (r)/dr  evaluated  at  r — 0 and 


exp  (//Co/4  (w/ (r) ) 1, 


A(w^(r))  -y  / ‘^''11“^  - /3^a^[w^*(r)  - p 


The  path  w/fr)  satisfies  the  differential  equation 


-p^aHwl  - P^^'h^) 
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with  iv.*(0)  = /3  and  w,*  (/?)  = /3  The  solution  to  Eq.  (8.47)  with  w/(0) 

fi  and  pCO)  =dw,{0)/dris 


Hence 


=P  + (fs  - ^A/)  cos  a/Jr)  + sin  aftr. 

{R)  =^~'^^[z^  + (Zj  —Z/^)cosapR]  + sin  apR. 


From  this  last  equation  we  have 


dyv,'(R)  sin  aPR 
dp  (0  ) ap 


The  caustics  are  at  apR  = nn,  n =1.2,  ....  If  Eq.  (8.49)  is  now  used  to  eliminate  p(0) 
from  the  expression  for  w,  (r),  it  only  takes  some  algebra  to  show 

A(w'(r))  = -~ap  H(z  -Za/)^  + (Zj  - z^)^]  cos  apR 
“ 2 Sin  fxpR 


- 2 (z  - z^)  (Zj  - z^)l. 


Collecting  terms  gives 


, i/z  .. 

dJ,  = , exp  -Z;w)^lcos«^7? 

^ 2iris\napR  2 sm  a/3/< 


-2(z  -z^)(z,  “^W>j|' 

The  solution  to  the  parabolic  equation  for  the  parabolic  profile  is  thus  extremely  simple. 
Although  this  solution  has  been  known  for  at  least  30  years,  it  is  continuously  being 
rediscovered.  Equation  (8.52)  may  be  generalized  to  include  the  effects  of  a pressure-release 
surface  by  using  the  techniques  of  Ref.  15,  and  a solution  may  also  be  readily  obtained  if  a has 
a range  dependence. 

Our  expression  for  4>j  may  be  written  as  an  infinite  sum  of  products  of  Hermite  polyno- 
mials //„  by  using  the  formula  179] 

" p-'""  f 1 ( iO  + 

exp]— ^ 


■f  „ ! „ [(x^  -1-  y^)  cos  0 - 2xyl  . 
2 sm  0 


Letting  9 — apR  gives 


<t>,  -p^''^  ZZ„(z  -Zm)Z„(z,  -zm)c^P 


—iaPR 


{2n  + \ ) . (8.53) 


where 


Z„  (z)  - |-^^|  (2  V ) (Vv7z)exp  | - | *„«z2j. 
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These  functions  obey  the  orthonormality  condition 

oo 

ZJz)Z,„(z) 

We  substitute  Eq.  (8.53)  into  Eq.  (8.41): 


p(x)  = ^ Z„  (z  - )Z„  (z,  - ) -^  / -^  exp 

n-O  0 P 


where 


ik„R 


1 

^ -f 

li  kl 


^ ~ ‘“^o  ^2”  + 1 )■ 

By  letting  fi  = k„T/k„,  the  integral  in  Eq.  (8.56)  becomes 


2 


1 


T + — 
T 


1 r dr 
— I — exp 

47r  -fi  T 

which  according  to  Eq.  (A  12)  is  (;/4)//J'’  {k„R).  Hence  we  obtain 

=T  -'^m)Z„(z,  -z^)H.l^Uk„R). 


(8.55) 


(8.56) 


(8.57) 


(8.58) 


which  is  the  standard  result. 


8.4.  Modified  Perturbation  Theory  and  the  Supereikonai  Approximation 

We  shall  consider  now  the  application  of  ['radkin’s  modified  perturbation  theory  (12)  to 
the  iiclmholtz  equation.  This  scheme  has  been  developed  for  underwater  sound  by  Munk  and 
Zachariasen  [80]  and  by  Callan  and  Zachariasen  [811  by  working  order  by  order  in  standard 
perturbation  theory.  Modified  perturbation  theory  gives  in  lowest  order  an  approximation, 
coined  the  "supereikonai  approximation"  in  Refs.  80  and  81,  closely  related  to  the  Rytov  ap- 
proximation [74-76]. 


Returning  to  Eqs.  (8.2)  and  (8.3),  we  define  ^ by 


where 


<1>  (t,  x|0,  X, ) “ <(>„  (t,  x|0,  ) exp  Q. 


d'„(T,  x|t  '.  x') 


K 

3/2 

ik„  (t  — r ') 

X - x' 

2 

Irri  (t  — t ') 

C X p 

I 

T — T ' 

tB.5'^) 


(8.60) 


The  function  d*,,  is  the  solution  to  Eqs.  (8.3)  for  « ~ I ( d/  - 0).  Equation  (8.3b)  will  be 
satisfied  if  y "*  0 at  r — 0.  It  follows  from  Eq.  (8.3a)  that  y satisfies 

- ^'k„H,Q — (x  - X,)  ■ vy  -t-  V^y  + (vy)2  - 2kjj  (/(%).  (8.61) 


We  wish  to  obtain  a series  expansion  for  yin  powers  of  the  perturbation  U(x).  The  simplest 
way  of  obtaining  this  is  to  replace  d/(x)  by  gl/ix),  where  p is  an  expansion  parameter  which  is 
set  equal  to  unity  at  the  end  of  the  calculation.  We  then  write 


y - i 

/-I 


(8.62) 
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After  substituting  Eq.  (8.62)  into  Eq.  (8.61),  we  equate  like  powers  of  g and  obtain  the 
coefficients  The  series  in  Eq.  (8.62)  starts  with  / = 1 because  <l>  = =0)  if 
fy  = 0.  For  / •=  1 we  have 


-2/MrP''’ 


(*  - X,)  • Ve'*'  + - 2k^U. 


The  terms  with  / > 1 are  obtained  from  the  recursion  relation 

2ik 

^ (x  -X,)  • 

^ ' m -I 

If  we  could  solve  Eqs.  (8.63)  and  (8.64),  we  would  have  the  exact  solution 

oe 

<J>  - '1‘oexp  £(?">.  (8.6S) 

/-I 

Here  we  will  be  content  with  the  approximation 

4>  - <J>„exp(?‘'>.  (8.66) 

The  higher-order  terms  have  been  discussed  by  Fradkin  (12). 

Equation  (8.63)  is  not  difficult  to  solve.  One  writes 

(?'  ” xtO,  X,)  (8.67) 

and  then  solves  the  resulting  differential  equation  for  Q'  by  using  standard  Green's-function 
techniques.  The  result  is 

° (t^x^O.  X, ) J I " '■  ’‘'^***» 

Comparison  of  Eqs.  (8.66)  and  (8.68)  with  Eq.  (6.81)  shows  the  approximation  Q =0*”  is 
equivalent  to  the  Rytov  approximation  for  dJ.  The  preceeding  analysis  represents  a different 
technique  for  obtaining  the  Rytov  approximation;  one  based  on  the  differential  equation 
satisfied  by  <t>  rather  than  the  path  integral. 

The  expression  for  may  be  cast  into  a different  form  by  introducing  the  Fourier 
transform 

U(x)  = f ’‘t/(k)  (8.69) 

(27r)^ 

and  carrying  out  the  integration  over  x'in  Eq.  (8.68); 

qII  ) « _ ii^  r — ) r dr 'exp  ik  ■ x + ~ (x  — x ) 

(27r)^  *5  ' T 

|l  ““1 

2k„  T ) 

(This  equation  is  to  be  compared  with  Eq.  (6.75b).)  Substituting  Eqs.  (8.60),  (8.66),  and  (8.70) 
into  Eq.  (8.2)  gives  the  supereikonal  expression  for  the  pressure  (80,81): 


1 

1/2 

f dr 

. 1* 

4n 

2nf 

2 

T “T 

r 

- 2t/(t)  + /*T  . (8.71) 
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where 


liT)  = f 


dk 


(lir ) ■ 


I 

t/(k ) J”  dli  exp 


/k  Ix^  +/a(x  -xj]  -/i)}.(8.72) 


The  Rylov  approximation  for  <1>  led  to  the  supereikonal  approximation  for  p(x).  On  the 
other  hand,  the  Rytov  approximation  applied  directly  to  p(\)  gives 


p(x)  =G„(x  -X,  )exp 


-2k 


G’„(x 


— r r dx' G’(x  — x ) (y(x')G'(x'  — X.) 

— X, ) ' 


(8.73) 


with 


G’  (x) 


47r  X 


i ^,'*,,1*1 


(8.74) 


(The  exponent  in  Eq.  (8.73)  is  the  Born  term  divided  by  G„.)  We  now  ask:  What  is  the  rela- 
tionship between  the  Rytov  approximation  for  pix)  and  the  supereikonal  approximation?  To 
answer  this  question,  we  return  to  Eq.  (8.68)  for  and  observe  that  the  integral  over  t' 

may  be  evaluated  with  the  aid  of  the  elementary  formula 


; 


dr ' 

lk„ 

I 

> 

. B^ 

(tTt“-  tTI^ 

2 

T — T ' 

T ■" — ~ 

T 

2nf 

kr 


1/2 


iAl_±_[Bi 

t|A||B| 


exp 


ik 


" (|A|  + |B|)^ 


(8.75) 


The  result  is 


A-rr 


4n- 
X exp 


lx  - x'  -t-  x'  - X, 


X - X X - X, 


~ [(|x  - x’l  -I-  |x'  - xj  )^  - lx  - x,l  ^1  . 


(8.76) 


In  the  straight-line  geometric-optics  approximation  the  integration  over  x'  is  restricted  to  the 
straight  line  joining  the  points  x^  and  x,  giving  |x  - x'|  -I-  |x'  - xJ  = |x  - xJ.  Here  we 
shall  assume  that  large-angle  scattering  is  unimportant,  so  that  the  difference  |x  - x'l  + 
|x'  — xJ  — |x  — x,|  is  small  compared  to  |x  — xJ.  With  this  assumption  Eq.  (8.76)  be- 
comes 


C>") ^ / ‘At 


|x  -xJ 


X - X X - X, 


X exp 


/A„|x  -X, 


(|x  - x1  -I-  |x  - xJ  - |x  - xJ  ) 


(8.77) 


Let  us  return  to  the  expression  for  p(x)  in  the  supereikonal  approximation: 


P(x)  - 


1 

i. 

1/2  « 

"o 

C dr 

4n 

Ini 

exp 


2 


r -f- 


|x  -xJ' 


+ Q 


(I) 


(8.78) 
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The  integral  over  t is  now  evaluated  by  the  method  of  stationary-phase  (Appendix  B).  In 
determining  the  stationary  phase  point,  the  dependence  of  Q*'*  on  t is  ignored.  We  find 


pfx) 


-*,l 

47t|  X — X J 


exp  'U 


(8,79) 


By  using  Eg.  (8.77)  for  we  readily  obtain  the  Rylov  approximation,  Eq.  (8.73).  One  may 
also  obtain  Eq.  (8.73)  directly  from  Eqs.  (8.76)  and  (8.78)  by  using  a saddle-point  approxima- 
tion to  evalute  the  integral  over  t. 


We  conclude  this  section  by  noting  the  straight-line  geometric-optics  approximation  for 
p(x)  is  easily  obtained  from  the  supereikonal  approximation.  We  drop  the  term  in  Eq. 
(8.72)  for  /(t)  and  find 


I 

/ •=  J*  d/3  t/(Xj  -E  /3  (x  - x, ) ) 


(8.80) 


independent  of  t.  Substituting  this  into  Eq.  (8.71)  and  evaluating  the  integral  over  t by  the 
method  of  stationary  phase,  we  get 

1»  - »,l  j 

p(x)  =G„(x  -Xj)exp  f dt  u\ 


X,  + I 


(x  - X, ) 


(8.81) 


In  obtaining  Eq.  (8.81),  we  have  considered  only  the  term 


in  determining  the  stationary-phase  point.  (Of  course  in  this  case  the  integral  over  t can  be 
evaluated  exactly.) 


9.  THE  PARABOLIC  APPROXIMATION 


9.1  Preliminary  Discussion 

We  have  derived  path-integral  representations  for  the  Helmholtz  equation  and  for  the 
two-  and  three-dimensional  parabolic  equations.  In  this  section  we  will  use  these  representa- 
tions to  discuss  the  parabolic  approximation. 


To  see  how  one  might  develop  the  parabolic  approximation,  let  us  consider  the  expres- 
sion for  p(x}  in  the  straight-line  geometric-optics  approximation,  Eq.  (8.81): 


p(x)  *=  f/„  (x  — x^ ) exp 


■«J 


-iA„  f dr  U x^  + r 


(x  - Xj ) 


X - X, 


(9.1) 


We  derived  this  expression  using  a two-step  process:  We  made  a straight-line  geometric-optics 
approximation  to  of  Eq.  (8.2);  then  we  evaluated  the  integral  over  r by  using  the  method  of 
stationary  phase,  determining  the  stationary-phase  point  from  that  part  of  the  phase  in  Eq.  (8.2) 
which  does  not  depend  on  U. 
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We  now  write  Eq.  (9.1)  in  the  form 

p(x)  “0’„(x  - *J  ). 


where 

ill  (r)  — exp 
Now  >li  (r)  satisfies  the  differential  equation 


f (x  - X J 

ik„  J dr  u X, 

II  I * ~ * 1 1 


where 


in  which 


— 2i/c„d,i/<  (r)  — U„|n'(r)  -l||i/i(r). 


nir)  - cjc  X,  4 (x  x,  )| 


fi)  = y/  (x  - ^ 4-  (>  4 ^ 

We  compare  these  equations  with  the  two-dimensional  parabolic  equation  (Ivq  (3  5)) 
-2(A„9,i;<  (r,  z)  “id/  4-A/|«^(r.z)  - I ||t/<  (r.  z ), 

where 

n(r.z)  ”c„Alr,  + ~ (x  -x^).y^  4-  (y  -yj.z. 

I "2  ^2 

in  which 

^2  = V ~ -*^1 ) ^ ~ 4' J 

and  with  the  parabolic  equation  in  three  dimensions  (Eq.  (3.12a)) 

-2ik^,d,i/>  (r.  y ,z)  ” Id/  4-  d/  4-  (r,  v.  z)  - I lli/<  (r.  z). 

where 


»(r.y,z}  - cJc 


-x,),>i.z, 


in  which 


Aft  = 


(9.2) 


(9.3) 


(9.4) 


(95a) 


(9.5b) 


(9.6) 


(9.7a) 


(9.7b) 


(9.8) 


(9.9a) 


(9.9b) 

The  .similarity  between  these  three  sets  of  equations  is  obvious  and  suggests  the  following 
scheme; 

• To  derive  a "parabolic  equation”  with  no  transverse  second  derivatives,  we  calculate  <l> 
using  straight-line  geometric  optics  and  evaluate  the  integral  over  r using  the  method  of  sta- 
tionary phase. 

• To  derive  a parabolic  equation  with  one  transverse  second  derivative,  we  calculate  4>  us- 
ing a straight-line  geometric-optics  approximation  applied  to  two  coordinates  and  then  evaluate 
the  integral  over  r as  before. 
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• To  obtain  a parabolic  equation  with  two  transverse  second  derivatives,  we  apply 
straight-line  geometric  optics  to  one  coordinate  and  then  make  the  stationary-phase  approxima- 
tion. 

In  the  next  subsection  we  will  demonstrate  the  correctness  of  this  scheme.  We  will  use  the 
first  path-integral  representation,  Eq.  (8.15),  of  the  previous  section.  Along  the  way  we  will  dis- 
cuss the  concept  of  additivity  of  the  action. 

A standard  reference  on  the  parabolic  approximation  is  the  paper  by  Klyatskin  and  Tartar- 
skii  [48].  In  subsection  9.3  we  will  review  their  calculation  and  compare  it  with  the  calculation 
in  subsection  9.2.  Klyatskin  and  Tatarskii  were  interested  in  the  application  of  the  parabolic 
equation  to  problems  in  the  theory  of  propagation  through  a random  medium.  To  discuss  this 
work,  we  must  consider  a model  for  the  random  fluctuations  on  the  sound  speed.  The  first  part 
of  subsection  9.3  will  be  devoted  to  this  model  and  to  an  approximation  intimately  related  to 
the  parabolic  approximation  called  the  Markov  approximation  [76,  Ch.  Sj.  In  subsection  9.4  we 
will  derive  an  improved  two-dimensional  parabolic  equation  by  relaxing  the  geometric-optics 
approximation.  In  the  final  subsection  we  will  discuss  corrections  to  the  stationary-phase 
approximation. 


9.2  Straight-Line  Geometric  Optics  and  the  Method  of  Stationary  Phase 


We  will  consider  first  the  derivation  of  the  two-dimensional  parabolic  equation.  According 
to  the  scheme  outlined  in  the  previous  subsection  this  requires  applying  straight-line  geometric 
optics  to  the  horizontal  coordinates.  Referring  to  Eq.  (9.11),  we  set 
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=.Vj^  + -y  +X^(fr) 


zUt)  = z'Ur).  (9.16) 

Here  we  are  using  a x subscript  to  indicate  a two-dimensional  vector  in  the  horizontal  plane: 

Ur  ) = (x(<r  ),  >'((r  ) ). 

R = ^ (x  - Xj ) ^ -t-  (j'  - }\ ) ^ = |xj^  -Si|.  (9.17) 

etc.  If  we  further  define  a two-dimensional  unit  vector  in  the  horizontal  direction  of  propaga- 
tion 


we  have 


ej,  ~^s-y  “7'j)  = (Jfi  - s^)/R, 


x^  Ur  ) = — r’i,  + Xj^  (<r). 


Equation  (9.11)  now  reads 

«I>(t,  x|0,  X ) = exp  — ^ f £)lx'(fr ) 1 exp  f d(T 

2t  •'  2 •'  dir 

X exp  —ikgj  dirU(s^  + x^^Ur),  z'Ur))  . (9.20) 

Horizontal  straight-line  geometric  optics  follows  from  ignoring  the  dependence  of  U on  x^Ur). 
With  this  approximation  the  integral  over  Xj^  (o-)  is  trivial,  and  we  obtain 

k \ ik 

<1>(t,  x|0,  X,)  *“  ■.  — exp  — ^ — (’(t,  zjO,  z ),  (9.21) 

l-trn  I 2t  I 


where 


r(T,  z|0,  Zj)  — J*  D(z(<r)lexp  ikg  ^ dir  yj^"!  — fi.  z(<r  )j  .(9.22) 


The  expression  for  the  pressure  is 


p(\)  “^^/“cxp  T + -^  r(T.rlO.  z ).  (9.23) 

4it  i)  T 2 [ T ) * 

We  now  can  learn  something  new  about  path  integrals.  One  might  expect  us  to  use  Eq. 
(9.22)  to  derive  a parabolic  equation  for  r(r,  z|0,  z, ).  In  fact  the  path  integral  in  Eq.  (9.22) 
does  not  in  general  satisfy  a simple  parabolic  equation,  and  the  reason  is  due  to  an  important 
concept  called  the  additivity  of  the  action. 
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To  illusiralc  this  concept,  consider  the  form  of  the  path  integral  which  solves  the  parabol- 
ic equation 

'P  (r,  z\r',  z')  = J” (v)  1 exp  (z  (.s) ) ],  (9.24) 

where 

r 

Aizis))  •=  J ds  L(s,  z{s)).  (9.25) 

f ' 

The  most  important  characteristic  of  this  path  integral  is  that  the  action  A is  additive: 

r"  r 

A(z(s))  = J*  ds  L (s,  Z|  (5) ) + J ds  L (i,  Zj  (5) ) 

r’  r" 

= /(|  (Z|  (s))  -E  (Z2  (.v)),  (9.26) 

where  Zj  (.s)  = z(.s)  for  r'  ^ .v  < r"and  Z2  (.v)  = z(s)  for  r"  ^ s ^ r.  The  component  actions 
/4|  and  A2  are  identical  in  form  to  A.  They  differ  only  in  the  limits  of  integration  and  in  the 
endpoint  conditions  satisfied  by  paths.  It  is  this  additivity  which  leads  to  the  composition  law, 
the  parabolic  equation,  and  the  marching  algorithms  for  numerically  evaluating  the  path 
integral.  If  the  action  is  not  additive,  one  has  none  of  these  things.  Generally  speaking,  path 
integrals  with  nonadditive  actions  are  very  difficult  to  handle.  The  application  of  Feynman’s 
theory  to  the  problem  of  propagation  through  a random  medium  has  not  been  particularly  fruit- 
ful, until  recently,  because  one  always  encountered  path  integrals  with  nonadditive  actions  (82). 
One  of  Dashen’s  most  important  contributions  (51  is  that  he  showed  these  nonadditive  actions 
become  additive  if  one  makes  the  Markov  approximation,  an  approximation  which,  as  we  shall 
see,  is  often  no  stronger  than  the  parabolic  approximation.  With  additive  actions  Dashen  was 
able  to  develop  a rather  complete  theory  of  propagation  through  random  media,  one  which 
seems  to  be  in  excellent  agreement  with  experiment  183). 

To  place  this  discussion  in  context,  we  observe  that  the  path  integral  of  Eq.  (9.22)  has  the 

form 

T 

J*  /2l  (z  (»r  ) 1 exp  J"  dtr  L(<r,  z(ir),  t)  . (9.27) 

The  action  is  not  additive,  because  r is  present  in  L.  If  the  sound  speed  is  a function  of  just 
depth,  1.  is  independent  of  t,  the  action  is  additive,  and  one  easily  obtains 

-2/A„;),I  (t.  z|0.  zJ  = -2klU(z)\V(T.z\Q.z,).  (9.28) 

with 

I (0,  z|0,  z,)  =fi(z  -z,).  (9.29) 

To  find  the  differential  equation  satisfied  by  I in  the  general  case,  we  return  to  the  discrete 
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where  Zjv  = ^ and  ir  — r.  A slraighlt'nrward  cakulaiion  gives 
I -2/A„d.  + 2Aity(x)ll  (t,  rlO.  rj 


2kf,R  r 'r  \ 

= i^l’('r)lexp  (A„  i/ir  y 


r dtT  ' cr  'e^  ■ V,  t/  Lj^  + z ((r  J 

i I ^ J 


2 T;ir  +--e,.z(a) 


where  = (3^^.  9^).  If  U is  independent  of  x^,  then  Vj^  f/  = 0 and  Eq.  (9.31)  reduces  to 
Eq.  (9.28). 

The  best  way  to  find  I'  is  not  by  solving  Eq.  (9.31)  but  by  introducing  a more  general 
path  integral: 


l\(r,  2|0,  r,)  = f D(z(rr))exp  /A„  f d(T  ^ — U +^^e.,r(rr)  .(9.32) 

'’n  2 d<r  T 


This  path  integral  is  more  complicated,  since  it  has  an  additional  dependence  on  the  variable  t. 
However 

r^(T,z|0,  z, ) = r (t,  z|0,  Zj ),  (9.33) 

and,  since  the  action  in  Eq.  (9.32)  is  additive,  it  is  not  difficult  to  show  ( 'obeys  a simple  para- 
bolic equation: 

— 2/A„3,I  ,.  (/,  zjO,  Zj ) = —2A, -I-  zj|  I\.  (/,  z|  0,  z, ),  (9.34a) 


r.,  (0,  z|0,  z,)  = 8 (z  —2,). 


(9.34b) 


Unlike  differential  equations  (9.28)  and  (9.31),  when  Eqs.  (9.34)  are  being  solved  t is  treated 
as  a parameter.  We  prefer  a slightly  different  form  of  Eqs.  (9.34)  obtained  by  letting 


(9.35a) 


9',(r,z|0.z,)  =r,  y .210,2,  . 


Equations  (9.23),  (9.34)  now  are 


1 f </t  '^0  , ^ 

I — exp  — T’  2 + H'.(/<,  z|0, 

•S  T 2 T 


(9.35b) 


(9.36a) 


;Y  3,'!',  (r,  zjO,  Zj)  - (dj  -(•  A^  [n^(r.  z)  - 1 11  M',  (r,  2|0.  z,),  (9.36b) 
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with  n{r,  z)  = nis^  + re^,  z)  and 

'J',  (0,  z|0,  z,)  = 8(z  - zj.  (9.36c) 

The  following  points  are  worth  discussing. 

• Equations  (9.36)  were  first  derived  in  Ref.  84.  In  that  paper  Feynman  path  integrals 
were  not  used.  Instead  the  pressure  field  p(\)  was  expanded  in  a power  series  in  U using  stan- 
dard perturbation  theory.  Each  term  in  the  series  was  then  approximated  by  applying  the 
straight-line  geometric-optics  approximation  (the  eikonal  approximation)  to  the  horizontal  coor- 
dinates, and  the  series  was  resummed.  The  calculation  made  no  assumptions  about  the  depth 
variation  of  p(\)  or  of  the  sound  speed.  Hence  Eqs.  (9.36)  are  valid  regardless  of  the  depth- 
dependent  boundary  conditions  imposed  on  the  pressure. 

• In  the  equation  satisfied  by  4'^  the  variable  of  integration  t enters  as  a parameter;  the 
equation  is  solved  as  if  it  were  constant.  If  t = Eqs.  (9.36b)  and  (9.36c)  reduce  to  the  pro- 
pagator defined  by  Eqs.  (3.9)  rather  than  the  field  t//  defined  by  Eq.  (3.8).  The  point  is  that  the 
depth-dependent  initial  field  h(z)  plays  no  role  when  one  considers  a point  source,  that  is, 
when  the  Helmholtz  equation  has  a spatial  8-function  on  the  right-hand  side. 


• In  deriving  these  equations,  we  have  assumed  only  that  straight-line  geometric  optics 
may  be  applied  to  the  horizontal  coordinates.  If  the  sound  speed  does  not  vary  with  the  hor- 
izontal coordinates,  this  type  of  approximation  is  exact.  Therefore  Eqs.  (9.36)  should  be  exact  if 
the  sound  speed  is  a function  simply  of  depth,  this  is  easily  demonstrated.  If  c(x)  = c(z), 
then 

=T  Z2„(z)Z„(z,)//,/'W/(„/?),  (9.37) 

^ n 

where 

k^n^(z) 
dz^ 

These  modal  functions  (which  we  assume  possess  only  a discrete  spectrum)  are  orthonormal, 


Z„(z) 


ZAz). 


(9.38) 


and  complete. 


/(7zZ„(z)Z„(z) 

0 

^ Z„(z)Z„(z')  =8(z  - z'). 


From  Eq.  (A12)  we  have  then 

1 f <A 


ik„R 

X +-L 

2 

X 

Letting  X = k„T/Rk„  gives  Eq.  (9.36a)  with 

(r,  zjO,  z, ) = ^ Z„(z)Z„(Zj)  exp 


nr 


lk.R  **" 


Using  Eqs.  (9.38)  and  (9.40),  it  is  easy  to  see  Eqs.  (9.36b)  and  (9.36c)  are  satisfied. 


(9.39) 


(9.40) 


(9.41) 


(9.42) 
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• We  have  emphasized  the  difficulty  of  writing  a differential  equation  for  r(T.z|0, 
when  the  sound  speed  depends  on  x and  >.  However  for  a range-independent  sound  speed  we 
immediately  obtained  Eqs.  (9.28)  and  (9.29).  Hence  for  this  case  Eqs.  (9.28),  and  (9.29)  must 
have  a relationship  to  Eqs.  (9.36b)  and  (9.36c).  To  exhibit  this  relationship,  we  note  that,  if 
-t-  z)  — n^(z),  the  form  of  Eq.  (9.36b)  dictates  9',  can  depend  on  ronly  through 
the  product  rr  (Eq.  (9.42)): 

9',.  (r,  z|0,  Zj ) = a function  of  rr. 

We  therefore  have  the  symmetry  relation 

/■d,9',  (r.  z|0,  z j = Td,9',  (r,  z|0,  z,).  (9.43) 

Substituting  this  into  Eq.  (9.36b),  we  get 

<r,  zlO,  z^)  = Id?  -f-/t?[n^(z)  - I 1)9'^  (r,  z|0,  z^ ).  (9,44) 

We  now  let  r = K: 

-2/A„f)/l',  («,  zlO,  zj  = (f)?  + k},\n'^(z)  - I 11  'P,  («,  zlO.  zj.  (9.45) 

Because  9',  {R.  z|l).  z, ) =1  (t,  z|0.  z,),  Eq.  (9.45)  is  equivalent  to  Eq.  (9.28).  In  the  general 
case,  'Pj  will  not  be  simply  a function  of  n because  of  the  additional  dependence  of  the  index 
of  refraction  on  r Therefore 

(/■('),  - TO,  )'P,  (r,  z|0,  z, ) 

— some  function  due  to  the  range  dependence  of  the  sound  speed;  (9.46) 

that  is,  the  symmetry  is  broken  by  the  range-dependent  part  of  the  sound  speed.  Since  the 
right-hand  side  of  Eq.  (9.46)  is  nonzero  (and  complicated),  it  is  not  in  general  possible  to  write 
/>(x)  as  an  integral  over  t of  an  exponential  limes  an  object  which  obeys  a simple  parabolic 
equation  with  a first  derivative  in  t.  (In  subsection  9.5  we  will  use  Eq.  (9.46)  as  the  starting 
point  of  a calculation  of  the  corrections  to  the  stationary-pha.se  approximation  due  to  the 
range-dependent  part  of  the  sound  speed.) 


To  complete  the  analysis,  the  integral  over  t in  Eq.  (9.36a)  is  approximated  by  the 
method  of  stationary  pha.se  (Appendix  B).  This  approximation  requires  that  ;*>  1 and 
that  'E,  be  a smoothly  varying  function  of  t in  the  region  of  the  stationary-phase  point 
1 — R.  We  easily  obtain  the  parabolic  approximation 


pix ) 


J. 

4jr 


2ni 


k„R 


1/2 


z|0,  z,) 


(9.47) 


where  'E  satisfies  the  parabolic  equation  given  by  Eq.  (3.9a)  with  the  fi-function  ini- 

tial condition  as  given  by  Eq.  (3.9b). 


We  now  consider  the  derivation  of  the  three-dimensional  parabolic  equation.  Application 
of  the  straight-line  geometric-optics  approximation  to  only  the  x coordinate  in  Eq.  (9.1 1 ) yields 


k 1 

1/2 

'll  (t.  xjo,  ) - 

exp 

2ntT  1 

ikg  (jf  — A ^ ^ 


I , <T.  p|0,  ), 


(9.48) 
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where 


l\(/.p|0,p,)  = /oipCfr ) ]exp 


ik 


of 


dp 

d(T 


-u 


Xj.  H (x  — Xj  ),p(<T  ) 


(9.49) 


Here  pUr)  = (yia),  zUt)),  so  pit)  = p = {y,  z)  and  p(o)  = p^  = (>>5,25).  Now  I" 
satisfies  the  equations 


-2klu 

X,  -E  — (x  - xj,  p 

T 

and 


Defining 


-2//(„a,i\  (/,  p|o,  Pj)  = 

r;(o.p|o,P,) 


r;(r,  p|0,  Pj),  (9.50a) 


t = 


(x  -Xj) 


(r,  p|0.  pj)  = 

we  have  the  system  of  equations 


(x  - Xj  ) 


. p|0,  p, 


(9.50b) 

(9.51a) 

(9.51b) 


1 

2ni 

f dr 

/k„ 

(x  - X ) ^ 1 

1 * 

47r 

K 

i r 

2 

T •+■ 

T 

(x  - Xj,  p|0,  Pj), (9.52a) 


-2ik 


H(x  -X  )]  ('•.pIO.p,)  = IV^  -2klU(x,  + r,p)]'i’^{r,p\0.p,),  (9.52b) 


and 


(0,  p|0,  Pj)  =6*^^(p-p^), 

The  stationary-phase  approximation  is  straight-forward  and  yields 


P(x)  ""’Vfx  -x,,p|0.p,), 

ZKq 


(9.52c) 

(9.53) 


where  ~ satisfies  Eqs.  (3.12). 


As  one  final  topic  we  consider  the  significance  of  the  parabolic  approximation  for  propa- 
gation in  an  infinite,  homogeneous  ocean  where  U = 0.  We  have 


'K/?.  2|0,  2j)  » 

Substituting  this  into  Eq.  (9.47)  we  get 

1 


'‘0 

1/2 

ikoR 

2 

2iTiR 

exp 

2 

R 

(9.54) 


p(x) 


4nR 


exp  j/)c 


R + 


(z  -z^)' 


2R 


(9.55) 


This  equation  is  to  be  compared  with  the  exact  result 


p(x) 


1 


Att^/R^  -I-  (2  -2,)2 


exp  ikg-jR  ^ + (2  - Zj ) ^ (9.56) 
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Obviously  the  parabolic  approximation  requires 

2 

« 1 (9.57a) 

and 

4 

« 277.  (9.57b) 

The  results  for  the  three-dimensional  equation  are  similar. 

9.3  The  Calculation  of  Klyalskin  and  Tatarskii:  Random  Fluctuations  in  the  Sound  Speed 

We  consider  in  this  subsection  the  situation  in  which  temperature  fluctuations  produce  a 
random  component  in  the  sound  speed.  The  function  U is  then  a random  variable  which  may 
be  written  in  the  form 

U(x)  = < UU)  > + (Uix)  - < U(x)  > ),  (9.58) 

where  < . . . > represents  an  average  over  a large  number  of  measurements  of  U.  This 

decomposition  defines  deterministic  and  random  components 

Uj(x)  = < U(x)>  (9.59a) 

and 

U^(x)  = U(x)  - <U{x)>.  (9.59b) 

We  shall  use  a particularly  simple  model  for  U,/  and  U^.  The  ocean  is  far  too  complicated 
for  this  model  to  be  even  approximately  valid,  and  calculations  based  on  it  can  at  best  give 
only  the  correct  orders  of  magnitude.  We  consider  it  not  because  we  fell  it  will  lead  to  a 
description  of  acoustic  fluctuations  but  simply  so  we  may  connect  our  analysis  of  the  parabolic 
approximation  with  a previous  study. 

We  assume: 

• The  deterministic  component  is  a function  of  all  three  spatial  coordinates: 

Uj  - Uj(x.y,z).  (9.60) 

• The  random  component  is  a Gaussian  random  variable.  Hence,  all  moments  of  can 
be  expressed  in  terms  of  the  correlation  function 

fl(x;  x')  - < U,(x)U,(x')  > . (9.61) 

Of  course  has  a zero  mean: 

<f/,(x)>  -0.  (9.62) 

• The  temperature  fluctuations  are  homogeneous  but  not  necessarily  isotropic: 

B(x,x')  - fl(|jr  - x'l,  |>-  ->■'1,  U - z'l  ).  (9.63) 


• The  fluctuations  along  the  coordinate  directions  are  characterized  by  scale  paramenters  « 

L^,  L^.  That  is,  for  |x|  small,  i 
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1 

I 


I 


I 

L 


B(x:0)  =<rl 


, 1 

X 

' 1 

y 

' _ 1 

z 

2 

' -y 

2 

h 

2 

L, 

+ ... 

(9.64) 


where 

= fl(0;  0)  = < ((/,(x)2)  > . (9.65) 

• the  correlation  function  possesses  separability  along  the  direction  of  propagation.  That 
is,  if  propagation  is  along  the  x axis, 


S(lxl,  |>’|,  I2I  ) = B(lxl,0,0.  )fi(0,  |>-1,  |z|  ). 
It  will  be  convenient  to  rewrite  this  equation  as 

Bdxl.Ij'Mzl)  =UU|)B(|y|.Ul) 

where 

+ 00 

C(|x|)  =fl(|x|.0,0)  / / dxB(\x\.Q.O) 

—00 

and 

00 

BdyMzD  = fi(o.  Ij'I.  lz|  ) J <)x5(|x|,0.0) 


(9.66) 

(9.67) 

(9.68a) 

(9.68b) 


• From  dimensional  considerations  we  have 


B(O.O)  = (9.69) 

and 

00 

/ dx  x'^%{\x\)  ^ Ll”  (9.70) 

— 00 

This  last  relation  will  be  used  only  for  n small,  (n  = 1 or  2). 

When  one  has  a random  sound  speed,  the  parabolic  approximation  is  intimately  related  to 
another  approximation  called  the  Markov  approximation  [76,  Ch.  5).  The  Markov  approxima- 
tion states  that  in  the  calculation  of  the  average  of  any  physical  quantity  the  function  I,  of  Eq. 
(9.67)  may  be  replaced  by  a 8 function 


4 (|x  - x'l ) — 8 (x  - x'). 


(9.71) 


Before  discussing  conditions  under  which  the  Markov  approximation  is  valid,  let  us  consider 
what  it  means.  If  the  correlation  function  has  a Gaussian  shape,  then 


«(UI ) 


1 

_ 1 

21 

X 

2 

(9.72) 


This  expression  approaches  a 8 function  as  — 0.  Consequently  the  Markov  approximation 
in  a problem  is  a statement  about  the  size  of  in  relation  to  the  other  dimensions  which 
characterize  the  effect  of  fluctuations  in  the  pressure  due  to  fluctuations  in  the  sound  speed 
along  the  direction  of  propagation.  If  L„  is  small  compared  to  these  dimensions,  it  may  be  set 
equal  to  zero,  giving  Eq.  (9.71). 


In  the  theory  of  propagation  through  a random  medium,  the  conditions  under  which 
some  approximation  is  valid  are  determined  by  comparing  expressions  for  the  average  value  of 
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( 


I 


products  of  the  pressure  field.  That  is,  expressions  for  < p>  , < pp*>  , < pp>  , < ppp*> 

before  and  after  making  the  approximation  are  compared.  In  principle  this  comparison  should 
be  done  for  all  average  values  of  experimental  interest.  In  practice  it  is  usually  done  only  for 
the  average  value  of  the  field  <p>  and  for  the  mutual  coherence  function  <pp'>.  In  this 
report  we  shall  consider  only  <p>.  Consequently  we  will  not  be  able  to  derive  all  the  condi- 
tions associated  with  a given  approximation. 


The  standard  technique  |76,  Ch.  S]  for  determining  the  conditions  for  the  validity  of  the 
Markov  approximation  is  to  assume  first  that  the  three-dimensional  parabolic  approximation  is 
valid.  We  therefore  return  to  the  path  integral  which  solves  the  three-dimensional  parabolic 
equation 


^'(X.pIO.pJ  D{p(r)\ 


exp 


X 

1 

dp 

ik„  1 dr 

0 J 

0 

2 

dr 

- U(r.p(r)) 


(9.73) 


Here  we  assume  the  direction  of  propagation  is  along  the  positive  x axis,  X ■=  x — is  the 
propagation  distance,  p - Cv,  r),  andp^  - (>’j,z, ). 


We  shall  calculate  < 'I'  > by  assuming  the  Markov  approximation  and  then  by  assuming 
{ is  a sharply  peaked  function  of  its  argument  rather  than  a 6 function.  Comparison  of  the  two 
expressions  will  give  the  first-order  correction  to  the  Markov  approximation  and  hence  a condi- 
tion for  its  validity. 


Taking  the  average  of  Eq.  (9.73),  we  obtain 


< 'I'  > - J*  D(p(r)  1 exp 


X 

ik„  r dr 

1 

OJ 

0 

2 

dr 

where 


€ (p  (r) ) - < exp 


- U,i(r.  p(r)) 


A 

-ik„  J dr  Uf  (r,  p(r) ) 


«(p(r)).  (9.74) 


> . 


Since  is  a zero-mean  Gaussian  random  variable,  we  can  write 

. -k^ 

((f>(r))  -exp 


--A 


with 


X X 


A - jf  dr  J*  dr' fl(|  r - r'l . |>'(r)  - ylr')] , | z (r)  - z(r')|  ). 


0 0 

We  now  use  Fq.  (9.67),  and  Eq.  (9.76b)  becomes 

X X 

A — J*  dr  { dr'(  (|r  — r'|  ) Fir,  r‘), 
n 0 

where  Fis  introduced  to  simplify  the  notation: 

F{r.r')  = B(\y(r)  - y (r')| . | z (r)  -z(r')|). 


(9.75) 


(9.76a) 


(9.76b) 


(9.77a) 


(9.77b) 
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There  is  a standard  technique  for  introducing  sum  and  difference  coordinates.  Since  the 
integrand  is  invariant  under  an  interchange  of  r and  r',  we  have 


A r 

\ ~ 2 J dr  J dr'  C,(\r  - r'\)  Fir.  r'>. 


If  we  let  r'  — r — s, 


A r 

A “ 2 f </r  f { (|  s|  ) Fir,  r — s) 

*0  o 

X X 

- 2 J*  </s  5 (|  s|  ) J*  Fir.  r — s). 


If  we  now  let  /■  — r + s/2,  we  have 


~ X-i/2 

A -2  J dsii\s\)  / di  F , + -^.t  -j-  . 
0 52  2 1 


We  have  replaced  the  upper  limit  on  the  integration  over  s by  since  we  always  assume 
X » L,.  and  4(|s()  = 0 for  s > L„.  Assuming  the  Markov  approximation,  we  have 

oo  X 

A - 2 f rfs8(s)  / rfr  Fit.  i). 

0 0 

oo 

But  Fit.  l)  — Bio,  0)  — <r  ?/  and  f rfs  8 (s)  — 4-,  so 

i)  2 

A ~ XL^trl-  (9.79) 

We  now  relax  the  8-function  assumption  and  assume  { (U|  ) is  only  sharply  peaked  about 
s “ 0.  We  then  can  examine  the  integral  over  / in  Eq.  (9.78)  in  a power  series  about  s — 0. 
We  shall  keep  terms  through  s^.  W'e  first  observe  that 


- B \s-^  + Ois^)\.\s-^  + Ois^)\ 
I Z at  at 


5(0.0)  1 


1 2 1 ^ ^ 1 2 1 * ^ 
2 ^ Ly  di  2^  Lz  dt 


+ O(s^).  (9.80) 


Hence 


I rf,  f - ij  - XL,,rl  • - 7 - y vj  + O(s^). 


where  y is  ihe  functional 


. a ^ ^ a a 'T  . ^ 

a,  ^ A 

LI  X i dr  n X i dr 


Substituting  Eq.  (9.81)  into  the  expression  for  A gives 


oo  I no 

2XL.,nl  /</s4(|.v|)  l - y J*  rfv  { (|sl  ) s^  . 
0 ( ^ ^ 0 
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The  first  integral  is  1/2  since  the  s! R term  may  be  dropped  for  X » 
is  L^l2  according  to  Eq.  (9.70).  Hence 


-1/2 


The  second  integral 


(9.84) 


This  expression  gives  the  first-order  correction  to  the  Markov  approximation.  By  comparing 
Eqs.  (9.79)  and  (9.84),  we  see  we  already  have  an  interesting  result.  The  Markov  approxima- 
tion is  valid  provided  y <s:  1 that  is,  provided 


and 


J_ 

X 


X 


\ 

f dr 

2 

0 

dr 

Lx 

X 

r 

dz 

2 

L, 

I dr 

«: 

Tr 

Lx 

(9.85a) 


(9.85b) 


The  Markov  approximation  implies  the  slopes  of  the  Feynman  paths,  averaged  over  the  range  of  propa- 
gation in  an  rms  sense,  must  he  small  [5].  In  the  ocean  where  ^ L^.  = JOO  L,,  Eq.  (9.85b)  is 
far  more  stringent  than  Eq.  (9.85a).  Returning  to  Eqs.  (9.76),  we  have 


- XL^k 


2,1-2 
o"  u 


«(p(r))  = exp 
Substituting  this  expression  into  the  path  integral  gives 

^ D\f>{r)]  exp 


-1/2 


^ y 


< > 


exp 


-i.,)fdt 

2 'k 


(I  - 


dt 


exp 


X 

ik„  f dr  {r,  p (r) ) 


1^ 

dt 


where 


and 


2 


rrl. 


(9.86) 


(9.87a) 


(9.87b) 


Equation  (9.86)  cannot  be  approximated  without  making  assumptions  about  the  deter- 
ministic sound  speed.  The  fact  that  the  Markov  approximation  is  intimately  related  to  the  pro- 
perties of  the  deterministic  sound  speed  follows  from  Eqs.  (9.85).  The  most  important  Feyn- 
man paths  will  be  the  deterministic  paths  calculated  using  ray  acoustics.  Therefore  the  magni- 
tudes of  the  left-hand  sides  of  Eqs.  (9.85)  are  measures  of  the  refractive  characteristics  of  the 
deterministic  sound  speed.  Rather  than  analyzing  Eq.  (9.86),  which,  in  fact,  is  not  well 
defined,  we  shall  adopt  a rather  naive  approach.  Since  — 0 in  the  Markov  approxima- 
tion, we  shall  take  as  cur  validity  conditions  1 and  Oj  is 
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ty 

2 

u 

2 

’ 

l<„Lj,ir  fj  •«  min 


For  = Ly  = L,  = L Ihis  reduces  to  the  standard  result 

k„ljrl  « 1. 

If  we  assume  Uj  — 0 and  the  Markov  approximation  is  valid,  Eq.  (9.86)  becomes 

I 21 


< S'  > 


IniX 


ik„  X 


In  this  case  Eq.  (9.53)  gives 
< p(x)  > 


1 


47r  X 

where  Lp  is  the  damping  distance. 


exp 


exp  /*„ 


1 


f‘s 


exp 


-XL^kl  irl 


X + — ~ 


exp  ( - XI Li)), 


kl  irl 


(9.88) 


(9.89) 


(9.90) 


(9.91) 


(9.92) 


From  Eq.  (9.88)  we  see  the  Markov  approximation  requires  that  Ld  be  large  in  comparison  to 
the  acoustic  wavelength. 


Having  discussed  the  Markov  approximation,  we  are  now  in  a position  to  discuss  the  pa- 
per by  Klyatskin  and  Tatarskii.  The  following  approach  is  adopted  in  Ref.  48;  First,  the  mean 
index  of  refraction  is  taken  to  be  a constant.  This  corresponds  to  setting  Uj(\)  equal  to  zero. 
Second,  a point  source  is  not  considered.  Instead,  it  is  assumed  that  a held  with  some  complex 
amplitude  distribution  i/„(p)  is  incident  on  the  x — 0 plane  and  that  the  random  perturbations 
in  the  sound  speed  exist  only  in  the  half-space  x > 0.  Third,  the  path-integral  representation 
given  by  Eq.  (8.20)  is  used.  Fourth,  the  fluctuations  in  the  index  of  refraction  are  homogene- 
ous and  isotropic.  Fifth,  although  it  is  not  necessary  to  introduce  a correlation  length,  validity 
conditions  are  recorded  in  terms  of  such  a quantity.  Sixth,  the  Markov  approximation  is  as- 
sumed. Seventh,  the  average  field  and  the  mutual  coherence  function  are  considered.  Only 
the  results  of  the  calculation  of  the  mutual  coherence  function  are  recorded. 


In  our  analysis  we  shall  assume:  the  mean  value  Uj  is  zero  (th*.  general  case  will  be  con- 
sidered in  the  following  subsections),  the  radiation  is  due  to  a point  source,  the  model  outlined 
above  describes  the  fluctuations,  and  the  Markov  approximation  is  valid.  Moreover,  we  will 
calculate  only  the  average  field,  and  we  will  use  the  path-integral  representation  we  listed  as 
Eqs.  (9.10)  through  (9.15). 


Returning  to  Eqs.  (9.10)  and  (9.1 1),  we  have 


< p(x)  > 


2 k„  D 

Just  as  before,  we  redefine  the  paths: 


I 


exp 


ik„T 


< <t>  (t,  x|0,  X ) > 


(9.93) 


x(cr  ) 


X,  -E  — (x  — X, ) -Ex'  (it  ), 
T ’ 


y(iT)  -^'((r). 
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z((r)  “ 


Ik  ^ 

<1>(t,  x|0,  Xj ) — exp  ( Dlx'  Ut)]  exp  ikg  J dtr 

2t  •'  *5 


1 ^ 
2 dtr 


— t/  |xj  + — A”  + x'  (o- ),  p (o-  )j  , 
where  again  X = x — Xj.  Taking  the  average  yields 

< x|0,  Xj)  > —exp  Z?lx((r)] 

'^o  C , dx  ^ 1 

X exp  I d<T  — — - — j 

2 -i  dv  2 


^ \d, 

dx 

2 r' 

^ diT 

in  which 


T T ^ 

A — f da-  J da'  ^ \ ~ X — — AT  + x(f7  ) — x(f7  ’)[  F(a,  a '). 

0 0 T T 

Fia , a')  B{\yia)  — y(a')\,\z(a)  — z((r')|). 


The  prime  on  the  path  has  been  dropped  in  writing  these  expressions.  We  need  only 
remember  that  now  we  have  x(0)  — x(t)  —0  Following  the  steps  which  led  to  Eq.  (9.78),  we 
get 

''-i/2  1.,  I 1 ( 1 1 1 1 


A — 2 f f rff  ( I 4 -1-  X r -1-  — X t 

* T 2 


F t ^ . (9.98) 

2 2 


If  we  were  to  apply  straight-line  geometric  optics  to  the  x coordinate,  x(r  -1-  5/2)—  x(/  —5/2) 
would  be  dropped.  Klyatskin  and  Tatarskii  go  one  step  further  and  calculate  A to  first  order  in 
this  difference 


^ i X . 1 5 -S 

U-jT^+xr-l-y  - X t - ~ 


With  the  Markov  approximation,  Eq.  (9.71),  this  last  equation  reads 
i - -y  8(5)  -f  -y  6 '(5) 


. s 

s 

r ^ 

— X 

^ 2 

Therefore 


A — Ai  -I-  A 3 -F 


(9.100) 


(9.101a) 
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where 


and 


-~f  dsSis)  / dtF  t +j-.t  -j- 

X i 2 2 


* ' T -V2 

^2  -^5  dsh-is)  / d, 
^ 0 7n 


X t + -f 

— X / — “ 

F / -h  V , f 

£ 

2 

2 

2 

2 

(9.101b) 


. (9.101c) 


The  function  A|  represents  the  result  of  applying  straight-line  geometric  optics  to  the  x coordi- 
nate, and  A2  is  the  first-order  correction  to  this  approximation.  It  is  not  difTicult  to  show 

,2 

(9.102a) 


A,  - — 


and 


A,  ral  \ dss  h (s)  \ dt  ^ 

^ X ^ i •'o  dt 


0. 


(9.102b) 


The  right-hand  of  this  last  equality  follows  because  x(0)  -jc(t)  - 0.  Therefore  the  first-order 
correction  to  straight-line  geometric  optics  vanishes  if  the  Markov  approximation  is  valid. 

We  now  drop  higher  order  corrections  to  straight-line  geometric  optics,  indicated  by  the 
dots  in  Eq.  (9.101a),  and  substitute  Eqs.  (9.102a)  and  (9.102b)  back  into  Eq.  (9.96); 


<<I>(t,x|0,  X,)  > 


expi-^^2  - 


^2-2  r _2 
*0  U 


or 


<<I>(t,  x|0,  Xj)  > 


lirir 


IX 

3/2 


J"  Z)[x(o-)l  exp 


^]d<r 

2 

2 •& 

rfo- 

. '*0  I I 7 T ^ 

«PlVl‘  -'.I  -UTr. 


In  writing  this  equation,  we  have  used 

|x  - Xjp  - X^  -t-  (p  - Pj)^ 
and  Eq.  (9.92),  the  definition  for  Lj).  Equation  (9.93)  becbmes 


An  integral  of  the  form 

I - 


1 

<^0 

1/2  « 

C dr 

1* 

4tr 

2ir/ 

\ t3/2 

2 

T T 

r 

1 

*0 

1/2  0. 

C dr 

'*0  . ^2 

4ir 

27r( 

1 2 T 

+ /(t) 


(9.103) 


(9.104) 


(9.105) 


(9.106) 


may  be  written  as 

QO  C« 

/ “-7^  / dcf  /(a)  f -77^  exp  { — ig<r  -I-  iAkg  -Jl  ■¥  Iq/kg).  (9.107) 

^irA  t27r ; 

By  using  the  expansion 

VI  -1  +...,  (9.108) 
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In  summary,  the  condition  for  the  validity  of  the  Markov  approximation  is  the  same  as  one  of 
the  conditions  for  the  validity  of  the  parabolic  approximation.  Moreover,  with  reference  to  Eqs. 
(9.85),  the  parabolic  approximation  implies  the  slopes  of  the  Feynman  paths  must  on  the  average  be 
small. 

9.4  Relaxing  the  Straight-Line  Geometric-Optics  Approximation:  A Modified  Parabolic  Equa- 
tion 

Having  discussed  the  calculation  by  Klyatskin  and  Tatarskii,  we  now  consider  corrections 
to  the  parabolic  approximation  in  the  general  case  in  which  the  sound  speed  is  an  arbitrary 
function  of  x.  in  this  subsection  the  straight-line  geometric-optics  approximation  will  be  re- 
laxed, and  in  the  next  subsection  we  will  discuss  corrections  to  the  stationary-phase  approxima- 
tion. Rather  than  decomposing  (/into  deterministic  and  random  components  and  considering 
the  average  value  of  the  field,  we  will  use  the  unaveraged  expressions  developed  in  subsections 
9.1  and  9.2.  Moreover  we  consider  only  the  two-dimensional  parabolic  equation. 
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Horizontal  straight-line  geometric  optics  amounted  to  assuming  that  in  Eq.  (9.20) 
(/  Si  + Ci  + Xi  (a-),  z ((r)|  — (/  Si  -I-  Ci,  z'  (a)j. 

T 1^1 

This  approximation  may  be  relaxed  in  several  ways.  We  could  take 
i/|si  -I-  f’i  -I-  Xi  (er),  z'(a-)j  — (/jsi  -I-  Ci,  z'  (o-)j 


(9.114) 


-I- Xi(o-)  - Vi  i/|si  -I- Ci.  z' (a)j  (9.115) 

and  evaluate  the  integral  over  xi(o-),  using  an  equation  similar  to  Eq.  (6.71).  We  would  end 
up  with  a two-dimensional  parabolic  equation  with  1/  replaced  by  an  expression  involving  in- 
tegrals over  the  path  parameter  rof  U and  U^.  Here  we  would  like  to  discuss  a different,  more 
general  approach.  Rather  than  making  a horizontal  straight-line  geometric-optics  approxima- 
tion, we  will  make  a horizontal  Rylov  approximation. 

Equation  (9.11)  may  be  written  in  the  form 


<I>  (t,  x|0,  Xj)  - J*  £) [z ((T ) 1 exp  rfer  j-^l  <I>j^  (t,  Xj^jO,  Sj^),  (9.116) 

where 

1,  t ^ j2 

'^0  ^ rfo-  y “ L/(x^(<t).  z(<r))  . (9.117) 

This  integral  is  a functional  of  the  path  z(ar).  It  can  be  evaluated  using  the  Rytov  approxima- 
tion (subsection  6.4): 


<l>i  (t,  Xj^IO.s^)  - (t.  xJO,  s^)  exp 


Born  term 


•froj.  (t,  sj  ’ 


(9.118) 


where 


T 

Born  term  “ — I da-  J*  dx^  <l>„^  (t,x^|o-,Xj^)  U{x^  z(a))  (o-,Xj^|0,Sj^ ) (9.119) 

(t,  x^l T x^ ) - J D\xja)]  expj-^  f da  j-^j 

f *0  1 \i  ~ \xi.  - x'i^Y 

- |2./(7^)  | H 2 W^\  ■ 


Substituting  Eqs.  (9.119)  and  (9.120)  into  Eq.  (9.118)  yields 


4,,(r,xJO.s,)  - — exp  -y-  exp 


ikoR^  r 

— exp  -ik„  J da 

2t  *5 


2 ff  /O’ (t  — o ) 


, '"o' 

dXi^  U{x^,z{a))  txp  - — r Sj^  -Xj^ 

2o(t— o) 


+ — (x^  - Sj^) 
T 


(9.121) 
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With  Eq.  (9.121),  Eq.  (9.1 16)  becomes 


x(0,  Xj  ) 


2irh 


exp 


ik„R^ 


r(T,  z|o,  Zj), 


(9.122) 


where 


r(T,r|0.z,)  -/  D(z(rr)J  exp 


*0  / 


2"^  - (/,;^(rr,  T,z(o-)) 


(9.123) 


with 


k..T 


T.  zUt))  ^ ^ f dx  U{x  z{<t)) 

7 TT  UT  It—  A 


IniiT  (t  — (t) 

i k„  T 


exp 


2 <r  (t  — tr) 


A’i  - x^  + — R 


(9.124) 


These  equations  are  to  be  compared  with  Eqs.  (9.21)  and  (9.22). 


Since  the  action  in  Eq.  (9.123)  is  not  additive,  we  generalize  (imbed)  just  as  before  by  in- 
troducing 


!',  (r,  z|0,  z, ) “ J*  D(z((r)lexp 


t 

ik.  1 d(r 

1 

dz 

"'6 

2 

dtr 

-U^hUt,  t,  z((r  ) ) 


. (9.125) 


We  have 


and 


where 


I\(t,z|0.  Zj)  “r(T,z|0,  Zj)  (9.126a) 

-2  ik„  B,  r,  (/,  z|0.zj  - (e^  - 2 *2  T,  z)l  f;  (/,  zlO,  z,),  (9.126b) 

(9.126c) 


r,  (0,  zjO,  Zj)  — 6 (z  “Zj ). 
In  analogy  with  Eqs.  (9.35)  we  define 


'P,  (r.  z|0,  z,)  -f; 
and  obtain  the  system  of  equations 


(9.127) 


i k„ 


T -E 


R‘ 


(/?,  z|0.  z,),  (9.128a) 


where 


-2  ik 


(r/R) 


- a, 'P,  (r,  z|0,  Zj)  - Id/  - Ik]  Uljj  {r,z)]'i^(r,z\Q,z,).  (9.128b) 
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in  which 


'P,  (0.  zjO.  Zj)  -8(z  -Zj), 


(9.128c) 


k„R(R/T) 


Kgl\\i\/Tf  f I ' kg  R iR/r ) ,1 


(9.128CI) 


If  the  integral  over  r in  Eq.  (9.128a)  is  evaluated  by  the  method  of  stationary  phase,  we  have 


4 TT  R 


(9.129a) 


-2i  kgd,^(r,z\0.z,)  -{d}  +klln^ff{r.z)  - 1 1) 'P  (r,  z|0,  z^),  (9.129b) 

’P(0,z|0,z,)  -6(z  -z,),  (9.129c) 

and 

~2nirlR-7)  ^2r(^-r)  ^ O ' • 

(9.1 29d) 

This  system  of  equations  was  first  derived  in  Ref.  84  using  different  techniques.  Equa- 
tion (9.129b)  is  identical  in  form  to  the  usual  two-dimensional  parabolic  equation  with  the  ex- 
ception that  the  index  of  refraction  is  replaced  by  an  effective  index  of  refraction  obtained  by 
integrating  over  the  horizontal  coordinates.  The  dominant  contributions  to  this  two- 
dimensional  integral  come  from  a cigar-shaped  region  about  the  geometric-optics  path.  The  ex- 
tent of  this  region  is  determined  by  the  horizontal  gradients  of  n^(x)  and  typically  is  of  the 
order  of  10  kms.  This  impoved  equation  is  not  particularly  useful  when  modeling  the  effects  of 
sound-speed  fluctuations  which  are  horizontally  homogeneous  and  isotropic.  It  does  offer  an 
advantage,  however,  if  one  is  interested  in  modeling  the  eiffects  of  large-scale  anomalies  such  as 
eddies  or  fronts,  because  it  allows  one  to  incorporate  into  a two-dimensional  parabolic  equation 
three-dimensional  variations  in  the  sound  speed. 

9.5  Corrections  to  the  Stationary-Phase  Approximation 

In  calculating  corrections  to  the  stationary-phase  approximation,  we  shall  assume  that  the 
horizontal  straight-line  geometric  optics  approximation  is  valid.  The  extension  to  the  modified 
parabolic  equation  developed  in  the  preceeding  subsection  is  straight  forward.  We  consider 
therefore  Eqs.  (9.36). 

Since  the  stationary-phase  approximation  amounts  to  the  assumption  that  is  slowly 
varying  about  the  stationary-phase  point  t “ R,  corrections  may  be  obtained  by  expanding  9',. 
about  T — /?: 


(R.zlO.z,)  - Lx  (t  -R) (R.zlO.z,), 

1-0  P- 


(9.130a) 
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where 

(«,  z|0,  z,)  - a/  z|0,  z,)|  , 

Substituting  Eq.  (9.130a)  into  Eq,  (9.36a)  gives 


with 


p(x) 


//(z) 


47r 


k„R 


,1,  a °°  dI 

^ I ^ /,(*„/?)  9r('>  (/(.Z|0.Z,), 

/-o  '■ 


z 

2ni 


1/2 


00 

'*  f — (f  — 1 )'  exp 

1 t 


y (/  + / -') 


(9.130b) 


(9.131) 


(9.132) 


The  calculation  of  corrections  to  the  stationary  phase  approximation  reduces  to  the  calculation 
of  the  integrals  //  and  the  derivatives  Mr 


We  first  indicate  how  /,  may  be  calculated  for  any  /.  By  the  change  of  variable  f — r it 
is  possible  to  write 


where  A/ is  the  coefficient  in  Eq.  (9.132): 

N - 


dt 

— exp 

-7  (/  + r') 

(r  - 1 )'  + 

1 

1 

i 

2 

1 

z 

Ini 


1/2 


(9.133a) 


(9.133b) 


For  any  / the  expression  in  the  second  pair  of  brackets  in  Eq.  (9  133)  may  be  written  as  an  /th 
order  polynominal  in  / + (l/t): 


U - D'  + 


-1 

I 


I K 


m —0 


t 


With  this  expression  Eq.  (9.133)  becomes 


//(-)  -y  I 


f ".'ar/fexp 


y (/  + /-') 


Equation  (A  12)  gives 


I (.,  1 "> 


/,(z)  I 4L 


m -0  I 


a I — //  (’ ) (2) 


Since  the  Hankel  function  obeys  Bessel’s  equation 

+-X  + ij  -0, 

dz^  z r/z  j ® 

it  is  always  possible  to  write 

//„<•  > (z)  - b„  // <'  > (z)  - c„  4 « ) (z), 


dz' 


~b„H„^^Uz)  +c„//,"'  (z). 


m "o 


(9.134) 


(9.135) 
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thus  expressing  all  the  higher  order  derivatives  in  terms  of  the  first  and  second  derivatives. 
The  expression  for  //  now  becomes 

/,(z)  ^ 1 (z)l-  (9136) 

^ m-ol  ' 


This  expression  gives  the  exact  value  of  the  integral.  The  determination  of  the  coefficients 
b„,  and  c„  is  straightforward  for  any  / and  m (but  tedius  for  large  / or  m).  As  an  illustration 
we  have,  for  / < 2, 


^0.2  /fq*  - -2  Al-0 

a}  aJ-  -2 

and 

f»0  - 1 f>,  - 0 />2 1 

Co  “ 0 c,  - -1  Cj  -1/z, 

giving 

/o(z)  »/  77  yv (z), 

(9.137a) 

7,  (z) in  N 1 + /-|-  7/”^  (z), 

(9.137b) 

and 

7,  (z)  - 2 / 77  yV  1 + /•  + - ^ 7/”>  (z). 

‘ Z UZ 

(9.137c) 

Since  z “ kgR  is  always  large,  it  is  usually  sufficient  to  approximate  the  Hankel  function  by  the 
first  term  of  its  asymptotic  series; 

,7r//«^(z)  = yv-'.  (9.138) 

Substituting  Eq.  (9.138)  into  Eqs.  (9.137)  and  carrying  out  the  indicated  differentiations  gives 
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► 


I 


(/?.  2|0,  z,)| , (R,  z|0,  z,). 


and  Eq.  (9.140)  becomes 


p(x) 

4v 


2iri 


. 1 1/2 


e'*'’'*'l'(/?.z|0,  Zj)  (1  + A), 


where  A is  an  estimate  of  the  error  associated  with  the  stationary-phase  approximation: 


A = 


iR 
2k  „ 


Bl'i’iR,  z|0,z,) 


'^fiR,  z|0,Zj) 


(9.141) 


(9.142) 


(9.143) 


This  is  as  far  as  one  can  go  without  specifying  the  sound  speed  and  using  a model  solution  to 
calculate  A via 


^(/?,  z|0,  Zj)  = ^ Z„  (z)  Z„  (Zj)  exp 


— (/t^  -k'^) 
2k,  ^ ^ 


(9.144) 


(This  equation  is  to  be  compared  with  Eq.  (9.42).)  For  a range-independent  sound  speed  the 
straight-line  geometric-optics  approximation  is  exact.  Therefore  the  stationary-phase  approxi- 
mation is  equivalent  to  the  parabolic  approximation,  and  the  question  of  the  validity  of  the  par- 
abolic approximation  reduces  to  the  question  of  the  size  of  A.  In  this  case  one  does  not  need  to 
numerically  solve  the  parabolic  equation  in  order  to  discuss  the  error.  All  that  is  required  is 
values  for  the  normal-mode  eigenfunctions  and  eigenvalues. 


If  the  sound  speed  possess  a range  dependence,  that  is,  if  it  depends  on  the  horizontal 
coordinates  Xj^  — (x,  y),  then  Eq.  (9.141)  is  no  longer  valid.  The  breakdown  of  Eq.  (9.141)  was 
suggested  when  we  wrote  Eq.  (9.46)  and  may  be  inferred  from  Eq.  (9.31).  We  could  calculate 
the  correction  to  the  stationary-phase  approximation  directly  from  Eq.  (9.140)  by  numerically 
solving  for  'i'^i.R,  z|0,  Zj)  for  a range  of  values  of  t about  t = R and  then  numerically  com- 
puting the  second  derivative.  For  a typical  range-dependent  sound  speed,  however,  the  most 
important  correction  to  the  stationary-phase  approximation  (not  the  parabolic  approximation)  is 
again  given  by  Eqs.  (9.142)  and  (9.143)  [85].  It  would  be  desirable  however  to  have  a criterion 
which  determines  when  Eqs.  (9.142)  and  (9.143)  are  likely  to  be  modified  as  the  result  of  the 
range  dependence  of  the  sound  speed. 


It  is  not  difficult  to  show  the  solution  to  Eq.  (9.36b)  satisfies  the  integral  equation 
'F,  (r,  z|0,  z,)  - ^lr  W)  (,,  ^|0,  z^)  dr^  f dzy 

X 9'/°^  (r.  z|r,,  z,  ) f/(/-,.  z,  ) (zj,  zJO,  Zj).  (9.145) 

Here  satisfies  Eq.  (9.36b)  with  f/  — — y (n^  — 1 ) set  equal  to  zero.  Since 

'*'(8’’  (r,  z|  r,.  z,  ) - Or.  z\  try  ,Zy),  (9.146) 

it  follows  that 

ir 

'^IR  z|0.  Zj)  - Or.  z|0,  Zj)  - ikg  ^ dry  J dzy 

X 'V^^Utr.zlry.Zy)  f/|-y,  z,  I (/•,  z,  1 0,  z, ).  (9.147) 
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Therefore 

'l',^(«,z|0,z,)  - ft '(//?,  z|0,Zj),  (9.148) 

where  ft  is  defined  by  the  equations 

z ft'(r,  zlO.Zj)  -0  (9.149) 

and 

ft'(0,z|0,  Zj)  -6(z  -z,).  (9.150) 

The  function  ft  depends  on  t in  two  distinct  ways.  First,  t scales  the  range.  This  depen- 
dence persists  even  when  the  sound  speed  is  independent  of  the  horizontal  coordinates. 
Second,  ft  has  a t dependence  which  results  solely  from  the  range  dependence  of  the  sound 
speed. 

In  analogy  with  Eq.  (9.131)  we  have 

P(x)  /"T  -D'exp  / -1--^  ft'(/?,z|0,  Zj).  (9.151) 

‘trr  i-Q  '■  0 ‘ fJJ 

If  we  take 

ft'(/?,z|0,  z,)  ==  ft '-'(/?,  z(0,Zj)  ='!'(/?,  z|0,z,),  (9.152) 

Eq.  (9.151)  reduces  to  Eq.  (9.140)  (for  kgR  large).  In  developing  a criterion  which  will  deter- 
mine when  Eqs.  (9.142)  and  (9.143)  are  valid,  it  is  sufficient  to  consider  only  the  / =•  0 term  in 
Eq.  (9.151): 

I °r  dt  1 

— exp  t n'(«.z|0,z,).  (9.153) 

(If  the  lowest  order  term  is  in  error,  there  is  little  hope  for  the  higher  order  terms.)  With  the 
stationary-phase  approximation,  Eq.  (9.153)  reduces  to  Eq.  (9.47)  because  of  Eq.  (9.152). 

We  now  write  a path  integral  for  ft: 

ft  '(/?,  z|0,  Zj)  ” J*  D[z(r)  \ exp  ikg  J*  ~ (9.154) 


The  t dependence  of  ft  occurs  only  in  the  argument  of  U.  Therefore  the  term  labeled  Eq. 
(9.153)  may  be  rewritten  in  the  form 
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By  carrying  out  an  analysis  similar  to  that  which  gave  Eq.  (9.110),  we  find 

1/2 


ylz(r)l  = 


1 


2iti 

k„R 


exp 


ik,,R  - ik„  ( dr  U(r,  zir))  exp  — - (f'U (r)  j ] 

0 


(9.157) 


where 


with 


R 

F[z(r)\  = k„  f drrU  (r.z(r)). 

■o 


(9.158a) 


U'(r.  z)  = 0,  U(r,  z). 


Substituting  Eq.  (9.157)  into  Eq.  (9.155)  yields 


1/2 

R 

Itri 

k„R 

^ D[z  (r)  1 exp 

ik„  ^ dr 

1 

2 

(Jz 

dr 

(9.158b) 


lUr.  z(r)) 


exp 


rill 

2k„R 


(9.159) 


Equation  (9.159)  agrees  with  Eq.  (9.47)  provided  /''^/2/c„/?  is  small,  say  a quarter  of  a cycle; 


2k„R 


< 


IT 

2 


or 


R 

f drr(/  (r,z(r)) 

*0 


< ttR. 


(9.160) 


(This  expression  is  analogous  to  Eq.  (9.113b).)  To  a good  approximation  Eq.  (9.160)  may  be 
replaced  by 


A 

/ drrU'(r,z\r)) 


< ttR 


(9.161) 


where  z (r)  is  any  one  of  the  ray  paths  defined  by  ray  acoustics.  Equation  (9.161)  gives  an 
operational  procedure  for  determining  when  the  leading  order  correction  to  the  stationary- 
phase  approximation  is  not  given  simply  by  Eq.  (9.143). 
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APPENDIX  A 

THE  FREE-SPACE  GREEN’S  FUNCTION  IN  n DIMENSIONS 

Let  3c  “ (x|, ....  x^)  be  an  n-dimensional  coordinate  vector,  and  let  “ X 

associated  Laplacian.  The  n-dimensional,  free-space  Green’s  function  G„  (xjx')  satisfies 

(V^  +/c^)C„(3c  -x') 6<''Ux  -3c'),  (Al) 

where  6 ^"^(3?)  is  the  «-dimensional  8-function.  In  this  appendix  a one-parameter  integral 
representation  for  G„  will  be  constructed. 


Introducing  the  Fourier  transform 

G„  (3c  - 3c')  - r e'*  -^'>G(Jt), 

" (2^.)" 


(A2) 


we  have 


cot) 


1 


— k\  — if 


(A3) 


for  outgoing  waves.  We  now  use 


oo 


</t  e 


-IT  (k}  —kl  —If) 


— kg  — if 

which  follows  from  the  Fourier  transform  of  the  step  function 

9(t)  “ — — r dut — , 

2ni  "oo  01  — If 

where 

e(T)  - 1,  T > 0, 

— 0,  T < 0. 


(A4) 


(A5) 


(A6) 


Substituting  Eqs.  (A3)  and  (A4)  into  Eq.  (A2)  and  interchanging  the  orders  of  integration  gives 


C.(J  - r)  - , I * /’*•' ; -0^  !'•  ''  • -'P' 


■) 


,ikj(x  —x')j  -ikjr 


(A7) 


/ I dr  e"""  n 

j-\  I 

The  n-dimensional  integral  over  k breaks  up  into  n one-dimensional  integrals,  and  each  of 
these  one-dimensional  integrals  can  be  readily  evaluated  (Eq.  (4.21))  giving 

1 n/2 

r Ilf 

G. 


7„(x  -X')  -i  j^dr 


irki  +iRllAj 


(A8) 
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where  R„  is  the  magnitude  of  x — x'; 


J-l 


(A9) 


With  the  change  of  variable 


2k„ 


we  have 


G„(x  -x')  - — 


1 

^0 

(n/2)  -1  » . 

r «T 

An 

2n/Rg 

4 ,W2 

2 

T 

T "T  — 

T 

(AlO) 


We  have  assumed  here  that  kg  is  real  so  that  k}  > 0.  Equation  (AlO)  remains  valid  if  kg  is 


complex  as  long  as  the  branch  cut  is  defined  so  that  kg  has  a positive  imaginary  part.  In  one 
dimension  with  “ |x  — x'|  we  have 


G,  (x  - x') 


/_ 

2k„ 


tk„R, 


1 27r/)?i  j J dT 


An  \ kg 


0 T 


1/2 


exp 


ikgRi  1 

— - — r + — 


(All) 


In  two  dimensions 


G,  (x  -X')  ~4-H}^HkgR,) 


oo 

1_ 

An  -I  T 


where 


0^2 

exp 


ikgR2 


T + - 
T 


(A12) 


and  in  three  dimensions 
Gj(x  - x')  ' 


y?2  ” "v/fx  — x')  ^ + (y 

1 

AnRj  " 

1 I kg  y'^7_dT 

['^3  ( ill 

An  [2ir;7l3  J 4 

1 2 r ^11 

(A13) 


with 


V (x  — x')^  + Cv  — y')^  + (z  —z')^. 


The  integral  representation  Eq.  (AlO)  may  be  used  to  derive  the  far-held  expression  for 
the  Green’s  function  G„.  If  kgR^  » 1,  the  integral  over  t in  Eq.  (AlO)  can  be  evaluated  by 
the  method  of  stationary  phase.  The  single,  stationary-phase  point  is  at  t — 1,  and  we  have, 
using  Eq.  (B4)  of  the  following  appendix. 


G„(x  -x')lrT-=r 


1 


AnR, 


2n  iR. 


€ 


(A  14) 
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THE  METHOD  OF  STATIONARY  PHASE 


Let 


l{p)  — J"  dr  ill 


(Bl) 


where  H is  a real  function  of  t with  two  continuous  derivatives  in  (a,  ft)  and  p is  some  nu- 
merically large  parameter.  If  0'(t,)  —0  (/  — 1, m)  with  a < t,  < /3  and  0"(t,)  ^ 0, 
then* 

' I 

ili{Ti)e’^  + 0\ 


Kp) 


I 

I -1 


2ir/ 


In  the  text  we  apply  the  method  of  stationary  phase  to  integrals  of  the  form 


(B2) 


/ - 

" 4 _n/2 


exp 


0 T 


K 

, b 

ar  + — 

2 

T 

>l>  (r). 


where  a and  bare  positive  parameters  and  is  assumed  large.  We  have 

(n  -3  )/4 


1 liri 

1/2 
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ik.yTab  , 
e ° ill 
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(B3) 
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